write_learning_fields.c

Go to the documentation of this file.
00001 /* ***************************************************** */
00002 /* Write learning fields for later use.                  */
00003 /* write_learning_fields.c                               */
00004 /* ***************************************************** */
00005 /* Author: Christian Page, CERFACS, Toulouse, France.    */
00006 /* ***************************************************** */
00011 /* LICENSE BEGIN
00012 
00013 Copyright Cerfacs (Christian Page) (2015)
00014 
00015 christian.page@cerfacs.fr
00016 
00017 This software is a computer program whose purpose is to downscale climate
00018 scenarios using a statistical methodology based on weather regimes.
00019 
00020 This software is governed by the CeCILL license under French law and
00021 abiding by the rules of distribution of free software. You can use, 
00022 modify and/ or redistribute the software under the terms of the CeCILL
00023 license as circulated by CEA, CNRS and INRIA at the following URL
00024 "http://www.cecill.info". 
00025 
00026 As a counterpart to the access to the source code and rights to copy,
00027 modify and redistribute granted by the license, users are provided only
00028 with a limited warranty and the software's author, the holder of the
00029 economic rights, and the successive licensors have only limited
00030 liability. 
00031 
00032 In this respect, the user's attention is drawn to the risks associated
00033 with loading, using, modifying and/or developing or reproducing the
00034 software by the user in light of its specific status of free software,
00035 that may mean that it is complicated to manipulate, and that also
00036 therefore means that it is reserved for developers and experienced
00037 professionals having in-depth computer knowledge. Users are therefore
00038 encouraged to load and test the software's suitability as regards their
00039 requirements in conditions enabling the security of their systems and/or 
00040 data to be ensured and, more generally, to use and operate it in the 
00041 same conditions as regards security. 
00042 
00043 The fact that you are presently reading this means that you have had
00044 knowledge of the CeCILL license and that you accept its terms.
00045 
00046 LICENSE END */
00047 
00048 
00049 
00050 
00051 
00052 
00053 
00054 #include <dsclim.h>
00055 
00057 int
00058 write_learning_fields(data_struct *data) {
00063   int istat; /* Diagnostic status */
00064 
00065   int ncoutid; /* NetCDF output file handle ID */
00066   int *timedimoutid; /* NetCDF time dimension output ID */
00067   int latdimoutid; /* NetCDF latitude dimension output ID */
00068   int londimoutid; /* NetCDF longitude dimension output ID */
00069   int sdimoutid; /* NetCDF season dimension output ID */
00070   int eofdimoutid; /* NetCDF EOF dimension output ID */
00071   int ptsdimoutid; /* NetCDF points dimension output ID */
00072   int *clustdimoutid; /* NetCDF clusters dimension output ID */
00073   int *weightdimoutid; /* NetCDF weight dimension output ID */
00074   int *timeoutid; /* NetCDF time variable ID */
00075   int latoutid; /* NetCDF latitude variable ID */
00076   int lonoutid; /* NetCDF longitude variable ID */
00077   int *cstoutid; /* NetCDF regression constant variable ID */
00078   int *regoutid; /* NetCDF regression coefficients variable ID */
00079   int *distoutid; /* NetCDF regression distances variable ID */
00080   int *rrdoutid; /* NetCDF precipitation index variable ID */
00081   int *rrooutid; /* NetCDF observed precipitation index variable ID */
00082   int *taoutid; /* NetCDF secondary large-scale field index variable ID */
00083   int *tadoutid; /* NetCDF secondary large-scale 2D field variable ID */
00084   int *rsqoutid; /* NetCDF regression R^2 variable ID */
00085   int *erroutid; /* NetCDF regression residuals variable ID */
00086   int *acoroutid; /* NetCDF regression autocorrelation variable ID */
00087   int *vifoutid; /* NetCDF regression VIF variable ID */
00088   int pcoutid; /* NetCDF pc_normalized_var variable ID */
00089   int tamoutid; /* NetCDF secondary large-scale field index mean variable ID */
00090   int tavoutid; /* NetCDF secondary large-scale field index variance variable ID */
00091   int *clustoutid; /* NetCDF clusters variable output ID */
00092   int *weightoutid; /* NetCDF weight variable ID */
00093   int vardimids[NC_MAX_VAR_DIMS]; /* NetCDF dimension IDs */
00094   
00095   size_t start[3]; /* Start element when writing */
00096   size_t count[3]; /* Count of elements to write */
00097 
00098   char *tmpstr = NULL; /* Temporary string */
00099 
00100   ut_system *unitSystem = NULL; /* Unit System (udunits) */
00101   ut_unit *dataunits = NULL; /* udunits variable */
00102 
00103   double fillvalue;
00104   float fillvaluef;
00105   char *nomvar = NULL;
00106   double *timeval = NULL;
00107   double *tancp_mean = NULL;
00108   double *tancp_var = NULL;
00109   double *bufd = NULL;
00110 
00111   int s;
00112   int t;
00113   int ii;
00114   int cov_true = FALSE; /* Check if cov is TRUE in at least one season */
00115 
00116   tancp_mean = (double *) malloc(data->conf->nseasons * sizeof(double));
00117   if (tancp_mean == NULL) alloc_error(__FILE__, __LINE__);
00118   tancp_var = (double *) malloc(data->conf->nseasons * sizeof(double));
00119   if (tancp_var == NULL) alloc_error(__FILE__, __LINE__);
00120 
00121   timedimoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00122   if (timedimoutid == NULL) alloc_error(__FILE__, __LINE__);
00123   clustdimoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00124   if (clustdimoutid == NULL) alloc_error(__FILE__, __LINE__);
00125   weightdimoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00126   if (weightdimoutid == NULL) alloc_error(__FILE__, __LINE__);
00127 
00128   timeoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00129   if (timeoutid == NULL) alloc_error(__FILE__, __LINE__);
00130   cstoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00131   if (cstoutid == NULL) alloc_error(__FILE__, __LINE__);
00132   regoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00133   if (regoutid == NULL) alloc_error(__FILE__, __LINE__);
00134   distoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00135   if (distoutid == NULL) alloc_error(__FILE__, __LINE__);
00136   rrdoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00137   if (rrdoutid == NULL) alloc_error(__FILE__, __LINE__);
00138   rrooutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00139   if (rrooutid == NULL) alloc_error(__FILE__, __LINE__);
00140   taoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00141   if (taoutid == NULL) alloc_error(__FILE__, __LINE__);
00142   tadoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00143   if (tadoutid == NULL) alloc_error(__FILE__, __LINE__);
00144   clustoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00145   if (clustoutid == NULL) alloc_error(__FILE__, __LINE__);
00146   weightoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00147   if (weightoutid == NULL) alloc_error(__FILE__, __LINE__);
00148   rsqoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00149   if (rsqoutid == NULL) alloc_error(__FILE__, __LINE__);
00150   erroutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00151   if (erroutid == NULL) alloc_error(__FILE__, __LINE__);
00152   acoroutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00153   if (acoroutid == NULL) alloc_error(__FILE__, __LINE__);
00154   vifoutid = (int *) malloc(data->conf->nseasons * sizeof(int));
00155   if (vifoutid == NULL) alloc_error(__FILE__, __LINE__);
00156 
00157   nomvar = (char *) malloc(200 * sizeof(char));
00158   if (nomvar == NULL) alloc_error(__FILE__, __LINE__);
00159   tmpstr = (char *) malloc(200 * sizeof(char));
00160   if (tmpstr == NULL) alloc_error(__FILE__, __LINE__);
00161 
00162   istat = utInit("");
00163 
00164   /* Open NetCDF file for writing, overwrite and truncate existing file if any */
00165   istat = nc_create(data->learning->filename_save_learn, NC_CLOBBER, &ncoutid);
00166   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00167 
00168   /* Set global attributes */
00169   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "processor", strlen(data->info->processor), data->info->processor);
00170   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "software", strlen(data->info->software), data->info->software);
00171   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "institution", strlen(data->info->institution), data->info->institution);
00172   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_email", strlen(data->info->creator_email), data->info->creator_email);
00173   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_url", strlen(data->info->creator_url), data->info->creator_url);
00174   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_name", strlen(data->info->creator_name), data->info->creator_name);
00175   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "contact_email", strlen(data->info->contact_email), data->info->contact_email);
00176   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "contact_name", strlen(data->info->contact_name), data->info->contact_name);
00177   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "other_contact_email", strlen(data->info->other_contact_email),
00178                           data->info->other_contact_email);
00179   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "other_contact_name", strlen(data->info->other_contact_name),
00180                           data->info->other_contact_name);
00181 
00182   fillvalue = -9999.9;
00183   fillvaluef = -9999.9;
00184 
00185   for (s=0; s<data->conf->nseasons; s++)
00186     if (data->conf->season[s].secondary_cov == TRUE) cov_true = TRUE;
00187 
00188   /* Set dimensions */
00189   istat = nc_def_dim(ncoutid, "season", data->conf->nseasons, &sdimoutid);
00190   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00191   istat = nc_def_dim(ncoutid, data->conf->eofname, data->learning->rea_neof, &eofdimoutid);
00192   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00193   istat = nc_def_dim(ncoutid, data->conf->ptsname, data->reg->npts, &ptsdimoutid);
00194   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00195   if (cov_true == TRUE) {
00196     istat = nc_def_dim(ncoutid, data->learning->sup_latname, data->learning->sup_nlat, &latdimoutid);
00197     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00198     istat = nc_def_dim(ncoutid, data->learning->sup_lonname, data->learning->sup_nlon, &londimoutid);
00199     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00200   }
00201 
00202   if (cov_true == TRUE) {
00203     /* Define lat and lon variables */
00204     vardimids[0] = latdimoutid;
00205     istat = nc_def_var(ncoutid, data->learning->sup_latname, NC_DOUBLE, 1, vardimids, &latoutid);
00206     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00207     istat = sprintf(tmpstr, "degrees_north");
00208     istat = nc_put_att_text(ncoutid, latoutid, "units", strlen(tmpstr), tmpstr);
00209     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00210     istat = sprintf(tmpstr, "latitude coordinate");
00211     istat = nc_put_att_text(ncoutid, latoutid, "long_name", strlen(tmpstr), tmpstr);
00212     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00213     istat = sprintf(tmpstr, "latitude");
00214     istat = nc_put_att_text(ncoutid, latoutid, "standard_name", strlen(tmpstr), tmpstr);
00215     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00216     
00217     vardimids[0] = londimoutid;
00218     istat = nc_def_var(ncoutid, data->learning->sup_lonname, NC_DOUBLE, 1, vardimids, &lonoutid);
00219     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00220     istat = sprintf(tmpstr, "degrees_east");
00221     istat = nc_put_att_text(ncoutid, lonoutid, "units", strlen(tmpstr), tmpstr);
00222     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00223     istat = sprintf(tmpstr, "longitude coordinate");
00224     istat = nc_put_att_text(ncoutid, lonoutid, "long_name", strlen(tmpstr), tmpstr);
00225     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00226     istat = sprintf(tmpstr, "longitude");
00227     istat = nc_put_att_text(ncoutid, lonoutid, "standard_name", strlen(tmpstr), tmpstr);
00228     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00229   }
00230 
00231   for (s=0; s<data->conf->nseasons; s++) {
00232 
00233     /* Define time dimensions and variables */
00234     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_time, s+1);
00235     istat = nc_def_dim(ncoutid, nomvar, data->learning->data[s].ntime, &(timedimoutid[s]));
00236     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00237 
00238     vardimids[0] = timedimoutid[s];
00239     istat = nc_def_var(ncoutid, nomvar, NC_INT, 1, vardimids, &(timeoutid[s]));
00240     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00241       
00242     istat = sprintf(tmpstr, "gregorian");
00243     istat = nc_put_att_text(ncoutid, timeoutid[s], "calendar", strlen(tmpstr), tmpstr);
00244     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00245     istat = sprintf(tmpstr, "%s", data->conf->time_units);
00246     istat = nc_put_att_text(ncoutid, timeoutid[s], "units", strlen(tmpstr), tmpstr);
00247     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00248     istat = sprintf(tmpstr, "time in %s", data->conf->time_units);
00249     istat = nc_put_att_text(ncoutid, timeoutid[s], "long_name", strlen(tmpstr), tmpstr);
00250     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00251 
00252     /* Define cluster dimensions */
00253     (void) sprintf(nomvar, "%s_%d", data->conf->clustname, s+1);
00254     istat = nc_def_dim(ncoutid, nomvar, data->conf->season[s].nclusters, &(clustdimoutid[s]));
00255     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00256 
00257     /* Define regression constant variables */
00258     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_reg_cst, s+1);
00259     vardimids[0] = ptsdimoutid;
00260     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &(cstoutid[s]));
00261     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00262 
00263     istat = nc_put_att_double(ncoutid, cstoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00264     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00265     istat = nc_put_att_text(ncoutid, cstoutid[s], "coordinates", strlen(data->conf->ptsname), data->conf->ptsname);
00266     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00267     istat = sprintf(tmpstr, "none");
00268     istat = nc_put_att_text(ncoutid, cstoutid[s], "units", strlen(tmpstr), tmpstr);
00269     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00270 
00271     /* Define regression coefficients variables */
00272     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_reg, s+1);
00273     vardimids[0] = clustdimoutid[s];
00274     vardimids[1] = ptsdimoutid;
00275     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 2, vardimids, &(regoutid[s]));
00276     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00277 
00278     istat = nc_put_att_double(ncoutid, regoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00279     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00280     istat = sprintf(tmpstr, "%s %s_%d", data->conf->ptsname, data->learning->nomvar_class_clusters, s+1);
00281     istat = nc_put_att_text(ncoutid, regoutid[s], "coordinates", strlen(tmpstr), tmpstr);
00282     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00283     istat = sprintf(tmpstr, "none");
00284     istat = nc_put_att_text(ncoutid, regoutid[s], "units", strlen(tmpstr), tmpstr);
00285     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00286 
00287     /* Define regression distances variables */
00288     if (data->learning->data[s].precip_reg_dist != NULL) {
00289       (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_reg_dist, s+1);
00290       vardimids[0] = timedimoutid[s];
00291       vardimids[1] = clustdimoutid[s];
00292       istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 2, vardimids, &(distoutid[s]));
00293       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00294       
00295       istat = nc_put_att_double(ncoutid, distoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00296       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00297       istat = sprintf(tmpstr, "%s_%d %s_%d", data->learning->nomvar_time, s+1, data->learning->nomvar_class_clusters, s+1);
00298       istat = nc_put_att_text(ncoutid, distoutid[s], "coordinates", strlen(tmpstr), tmpstr);
00299       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00300       istat = sprintf(tmpstr, "none");
00301       istat = nc_put_att_text(ncoutid, distoutid[s], "units", strlen(tmpstr), tmpstr);
00302       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00303     }
00304 
00305     /* Define regression R^2 diagnostic */
00306     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_reg_rsq, s+1);
00307     vardimids[0] = ptsdimoutid;
00308     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &(rsqoutid[s]));
00309     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00310     
00311     istat = nc_put_att_double(ncoutid, rsqoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00312     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00313     istat = nc_put_att_text(ncoutid, rsqoutid[s], "coordinates", strlen(data->conf->ptsname), data->conf->ptsname);
00314     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00315     istat = sprintf(tmpstr, "none");
00316     istat = nc_put_att_text(ncoutid, rsqoutid[s], "units", strlen(tmpstr), tmpstr);
00317     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00318     
00319     /* Define regression residuals diagnostic */
00320     if (data->learning->data[s].precip_reg_err != NULL) {
00321       (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_reg_err, s+1);
00322       vardimids[0] = timedimoutid[s];
00323       vardimids[1] = ptsdimoutid;
00324       istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 2, vardimids, &(erroutid[s]));
00325       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00326       
00327       istat = nc_put_att_double(ncoutid, erroutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00328       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00329       istat = sprintf(tmpstr, "%s %s_%d", data->conf->ptsname, data->learning->nomvar_time, s+1);
00330       istat = nc_put_att_text(ncoutid, erroutid[s], "coordinates", strlen(tmpstr), tmpstr);
00331       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00332       istat = sprintf(tmpstr, "none");
00333       istat = nc_put_att_text(ncoutid, erroutid[s], "units", strlen(tmpstr), tmpstr);
00334       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00335     }
00336 
00337     /* Define regression autocorrelation diagnostic */
00338     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_reg_acor, s+1);
00339     vardimids[0] = ptsdimoutid;
00340     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &(acoroutid[s]));
00341     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00342     
00343     istat = nc_put_att_double(ncoutid, acoroutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00344     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00345     istat = nc_put_att_text(ncoutid, acoroutid[s], "coordinates", strlen(data->conf->ptsname), data->conf->ptsname);
00346     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00347     istat = sprintf(tmpstr, "none");
00348     istat = nc_put_att_text(ncoutid, acoroutid[s], "units", strlen(tmpstr), tmpstr);
00349     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00350     
00351     /* Define regression VIF diagnostic */
00352     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_reg_vif, s+1);
00353     vardimids[0] = clustdimoutid[s];
00354     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &(vifoutid[s]));
00355     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00356     
00357     istat = nc_put_att_double(ncoutid, vifoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00358     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00359     istat = nc_put_att_text(ncoutid, vifoutid[s], "coordinates", strlen(data->learning->nomvar_class_clusters),
00360                             data->learning->nomvar_class_clusters);
00361     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00362     istat = sprintf(tmpstr, "none");
00363     istat = nc_put_att_text(ncoutid, vifoutid[s], "units", strlen(tmpstr), tmpstr);
00364     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00365 
00366     /* Define precipitation index variables */
00367     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_index, s+1);
00368     vardimids[0] = timedimoutid[s];
00369     vardimids[1] = ptsdimoutid;
00370     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 2, vardimids, &(rrdoutid[s]));
00371     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00372 
00373     istat = nc_put_att_double(ncoutid, rrdoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00374     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00375     istat = sprintf(tmpstr, "%s %s_%d", data->conf->ptsname, data->learning->nomvar_time, s+1);
00376     istat = nc_put_att_text(ncoutid, rrdoutid[s], "coordinates", strlen(tmpstr), tmpstr);
00377     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00378     istat = sprintf(tmpstr, "none");
00379     istat = nc_put_att_text(ncoutid, rrdoutid[s], "units", strlen(tmpstr), tmpstr);
00380     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00381 
00382     /* Define precipitation index obs variables */
00383     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_precip_index_obs, s+1);
00384     vardimids[0] = timedimoutid[s];
00385     vardimids[1] = ptsdimoutid;
00386     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 2, vardimids, &(rrooutid[s]));
00387     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00388 
00389     istat = nc_put_att_double(ncoutid, rrooutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00390     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00391     istat = sprintf(tmpstr, "%s %s_%d", data->conf->ptsname, data->learning->nomvar_time, s+1);
00392     istat = nc_put_att_text(ncoutid, rrooutid[s], "coordinates", strlen(tmpstr), tmpstr);
00393     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00394     istat = sprintf(tmpstr, "none");
00395     istat = nc_put_att_text(ncoutid, rrooutid[s], "units", strlen(tmpstr), tmpstr);
00396     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00397 
00398     /* Define sup_index (secondary large-scale field index for learning period) */
00399     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_sup_index, s+1);
00400     vardimids[0] = timedimoutid[s];
00401     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &(taoutid[s]));
00402     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00403 
00404     istat = nc_put_att_double(ncoutid, taoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00405     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00406     istat = sprintf(tmpstr, "%s_%d", data->learning->nomvar_time, s+1);
00407     istat = nc_put_att_text(ncoutid, taoutid[s], "coordinates", strlen(tmpstr), tmpstr);
00408     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00409     istat = sprintf(tmpstr, "none");
00410     istat = nc_put_att_text(ncoutid, taoutid[s], "units", strlen(tmpstr), tmpstr);
00411     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00412 
00413     /* Define sup_val (secondary large-scale 2D field for learning period) */
00414     if (data->conf->season[s].secondary_cov == TRUE && data->learning->data[s].sup_val != NULL) {
00415       (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_sup_val, s+1);
00416       vardimids[0] = timedimoutid[s];
00417       vardimids[1] = latdimoutid;
00418       vardimids[2] = londimoutid;
00419       istat = nc_def_var(ncoutid, nomvar, NC_FLOAT, 3, vardimids, &(tadoutid[s]));
00420       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00421       
00422       istat = nc_put_att_double(ncoutid, tadoutid[s], "missing_value", NC_FLOAT, 1, &fillvalue);
00423       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00424       istat = sprintf(tmpstr, "%s_%d %s %s", data->learning->nomvar_time, s+1, data->learning->sup_latname, data->learning->sup_lonname);
00425       istat = nc_put_att_text(ncoutid, tadoutid[s], "coordinates", strlen(tmpstr), tmpstr);
00426       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00427       istat = sprintf(tmpstr, "none");
00428       istat = nc_put_att_text(ncoutid, tadoutid[s], "units", strlen(tmpstr), tmpstr);
00429       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00430     }
00431   }
00432 
00433   /* Define pc_normalized_var */
00434   (void) strcpy(nomvar, data->learning->nomvar_pc_normalized_var);
00435   vardimids[0] = eofdimoutid;
00436   istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &pcoutid);
00437   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00438 
00439   istat = nc_put_att_double(ncoutid, pcoutid, "missing_value", NC_DOUBLE, 1, &fillvalue);
00440   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00441   istat = nc_put_att_text(ncoutid, pcoutid, "coordinates", strlen(data->conf->eofname), data->conf->eofname);
00442   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00443   istat = sprintf(tmpstr, "none");
00444   istat = nc_put_att_text(ncoutid, pcoutid, "units", strlen(tmpstr), tmpstr);
00445   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00446 
00447   /* Define tancp_mean */
00448   (void) strcpy(nomvar, data->learning->nomvar_sup_index_mean);
00449   vardimids[0] = sdimoutid;
00450   istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &tamoutid);
00451   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00452 
00453   istat = nc_put_att_double(ncoutid, tamoutid, "missing_value", NC_DOUBLE, 1, &fillvalue);
00454   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00455   istat = sprintf(tmpstr, "season");
00456   istat = nc_put_att_text(ncoutid, tamoutid, "coordinates", strlen(tmpstr), tmpstr);
00457   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00458   istat = sprintf(tmpstr, "none");
00459   istat = nc_put_att_text(ncoutid, tamoutid, "units", strlen(tmpstr), tmpstr);
00460   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00461 
00462   /* Define tancp_var */
00463   (void) strcpy(nomvar, data->learning->nomvar_sup_index_var);
00464   vardimids[0] = sdimoutid;
00465   istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &tavoutid);
00466   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00467 
00468   istat = nc_put_att_double(ncoutid, tavoutid, "missing_value", NC_DOUBLE, 1, &fillvalue);
00469   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00470   istat = sprintf(tmpstr, "season");
00471   istat = nc_put_att_text(ncoutid, tavoutid, "coordinates", strlen(tmpstr), tmpstr);
00472   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00473   istat = sprintf(tmpstr, "none");
00474   istat = nc_put_att_text(ncoutid, tavoutid, "units", strlen(tmpstr), tmpstr);
00475   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00476 
00477   /* End definition mode */
00478   istat = nc_enddef(ncoutid);
00479   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00480 
00481   /* Write variables */
00482 
00483   if (cov_true == TRUE) {
00484     /* Write lat and lon */
00485     start[0] = 0;
00486     start[1] = 0;
00487     start[2] = 0;
00488     count[0] = (size_t) data->learning->sup_nlat;
00489     count[1] = 0;
00490     count[2] = 0;
00491     istat = nc_put_vara_double(ncoutid, latoutid, start, count, data->learning->sup_lat);
00492     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00493     
00494     start[0] = 0;
00495     start[1] = 0;
00496     start[2] = 0;
00497     count[0] = (size_t) data->learning->sup_nlon;
00498     count[1] = 0;
00499     count[2] = 0;
00500     istat = nc_put_vara_double(ncoutid, lonoutid, start, count, data->learning->sup_lon);
00501     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00502   }
00503 
00504   /* Initialize udunits */
00505   ut_set_error_message_handler(ut_ignore);
00506   unitSystem = ut_read_xml(NULL);
00507   ut_set_error_message_handler(ut_write_to_stderr);
00508   dataunits = ut_parse(unitSystem, data->conf->time_units, UT_ASCII);
00509 
00510   timeval = NULL;
00511   for (s=0; s<data->conf->nseasons; s++) {
00512 
00513     timeval = (double *) realloc(timeval, data->learning->data[s].ntime * sizeof(double));
00514     if (timeval == NULL) alloc_error(__FILE__, __LINE__);
00515     
00516     /* Compute time variable using actual units */
00517     for (t=0; t<data->learning->data[s].ntime; t++)
00518       istat = utInvCalendar2(data->learning->data[s].time_s->year[t], data->learning->data[s].time_s->month[t],
00519                              data->learning->data[s].time_s->day[t], data->learning->data[s].time_s->hour[t],
00520                              data->learning->data[s].time_s->minutes[t], data->learning->data[s].time_s->seconds[t],
00521                              dataunits, &(timeval[t]));
00522     
00523     /* Write time */
00524     start[0] = 0;
00525     start[1] = 0;
00526     start[2] = 0;
00527     count[0] = (size_t) data->learning->data[s].ntime;
00528     count[1] = 0;
00529     count[2] = 0;
00530     istat = nc_put_vara_double(ncoutid, timeoutid[s], start, count, timeval);
00531     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00532 
00533     /* Write regression constants */
00534     start[0] = 0;
00535     start[1] = 0;
00536     start[2] = 0;
00537     count[0] = (size_t) data->reg->npts;
00538     count[1] = 0;
00539     count[2] = 0;
00540     istat = nc_put_vara_double(ncoutid, cstoutid[s], start, count, data->learning->data[s].precip_reg_cst);
00541     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00542 
00543     /* Write regressions coefficients */
00544     start[0] = 0;
00545     start[1] = 0;
00546     start[2] = 0;
00547     count[0] = (size_t) data->conf->season[s].nclusters;
00548     count[1] = (size_t) data->reg->npts;
00549     count[2] = 0;
00550     istat = nc_put_vara_double(ncoutid, regoutid[s], start, count, data->learning->data[s].precip_reg);
00551     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00552 
00553     /* Write regressions distances */
00554     if (data->learning->data[s].precip_reg_dist != NULL) {
00555       start[0] = 0;
00556       start[1] = 0;
00557       start[2] = 0;
00558       count[0] = (size_t) data->learning->data[s].ntime;
00559       count[1] = (size_t) data->conf->season[s].nclusters;
00560       count[2] = 0;
00561       istat = nc_put_vara_double(ncoutid, distoutid[s], start, count, data->learning->data[s].precip_reg_dist);
00562       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00563     }
00564 
00565     /* Write regression R^2 diagnostic */
00566     start[0] = 0;
00567     start[1] = 0;
00568     start[2] = 0;
00569     count[0] = (size_t) data->reg->npts;
00570     count[1] = 0;
00571     count[2] = 0;
00572     istat = nc_put_vara_double(ncoutid, rsqoutid[s], start, count, data->learning->data[s].precip_reg_rsq);
00573     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00574     
00575     if (data->learning->data[s].precip_reg_err != NULL) {
00576       /* Write regression residuals diagnostic */
00577       start[0] = 0;
00578       start[1] = 0;
00579       start[2] = 0;
00580       count[0] = (size_t) data->learning->data[s].ntime;
00581       count[1] = (size_t) data->reg->npts;
00582       count[2] = 0;
00583       istat = nc_put_vara_double(ncoutid, erroutid[s], start, count, data->learning->data[s].precip_reg_err);
00584       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00585     }
00586     
00587     /* Write regression autocorrelation diagnostic */
00588     start[0] = 0;
00589     start[1] = 0;
00590     start[2] = 0;
00591     count[0] = (size_t) data->reg->npts;
00592     count[1] = 0;
00593     count[2] = 0;
00594     istat = nc_put_vara_double(ncoutid, acoroutid[s], start, count, data->learning->data[s].precip_reg_autocor);
00595     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00596     
00597     /* Write regression VIF diagnostic */
00598     start[0] = 0;
00599     start[1] = 0;
00600     start[2] = 0;
00601     count[0] = (size_t) data->conf->season[s].nclusters;
00602     count[1] = 0;
00603     count[2] = 0;
00604     istat = nc_put_vara_double(ncoutid, vifoutid[s], start, count, data->learning->data[s].precip_reg_vif);
00605     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00606     
00607     /* Write precipitation index */
00608     start[0] = 0;
00609     start[1] = 0;
00610     start[2] = 0;
00611     count[0] = (size_t) data->learning->data[s].ntime;
00612     count[1] = (size_t) data->reg->npts;
00613     count[2] = 0;
00614     istat = nc_put_vara_double(ncoutid, rrdoutid[s], start, count, data->learning->data[s].precip_index);
00615     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00616 
00617     /* Write precipitation index obs */
00618     start[0] = 0;
00619     start[1] = 0;
00620     start[2] = 0;
00621     count[0] = (size_t) data->learning->data[s].ntime;
00622     count[1] = (size_t) data->reg->npts;
00623     count[2] = 0;
00624     istat = nc_put_vara_double(ncoutid, rrooutid[s], start, count, data->learning->data[s].precip_index_obs);
00625     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00626 
00627     /* Write secondary field index */
00628     start[0] = 0;
00629     start[1] = 0;
00630     start[2] = 0;
00631     count[0] = (size_t) data->learning->data[s].ntime;
00632     count[1] = 0;
00633     count[2] = 0;
00634     istat = nc_put_vara_double(ncoutid, taoutid[s], start, count, data->learning->data[s].sup_index);
00635     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00636 
00637     /* Write secondary field 2D field */
00638     if (data->conf->season[s].secondary_cov == TRUE && data->learning->data[s].sup_val != NULL) {
00639       start[0] = 0;
00640       start[1] = 0;
00641       start[2] = 0;
00642       count[0] = (size_t) data->learning->data[s].ntime;
00643       count[1] = (size_t) data->learning->sup_nlat;
00644       count[2] = (size_t) data->learning->sup_nlon;
00645       istat = nc_put_vara_double(ncoutid, tadoutid[s], start, count, data->learning->data[s].sup_val);
00646       if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00647     }
00648       
00649     tancp_mean[s] = data->learning->data[s].sup_index_mean;
00650     tancp_var[s] = data->learning->data[s].sup_index_var;
00651   }
00652 
00653   /* Write pc_normalized_var */
00654   start[0] = 0;
00655   start[1] = 0;
00656   start[2] = 0;
00657   count[0] = (size_t) data->learning->rea_neof;
00658   count[1] = 0;
00659   count[2] = 0;
00660   bufd = (double *) malloc(data->learning->rea_neof * sizeof(double));
00661   if (bufd == NULL) alloc_error(__FILE__, __LINE__);
00662   for (ii=0; ii<data->learning->rea_neof; ii++)
00663     bufd[ii] = sqrt(data->learning->pc_normalized_var[ii]);
00664   istat = nc_put_vara_double(ncoutid, pcoutid, start, count, bufd);
00665   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00666   (void) free(bufd);
00667 
00668   /* Write tancp_mean */
00669   start[0] = 0;
00670   start[1] = 0;
00671   start[2] = 0;
00672   count[0] = (size_t) data->conf->nseasons;
00673   count[1] = 0;
00674   count[2] = 0;
00675   istat = nc_put_vara_double(ncoutid, tamoutid, start, count, tancp_mean);
00676   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00677 
00678   /* Write tancp_var */
00679   start[0] = 0;
00680   start[1] = 0;
00681   start[2] = 0;
00682   count[0] = (size_t) data->conf->nseasons;
00683   count[1] = 0;
00684   count[2] = 0;
00685   istat = nc_put_vara_double(ncoutid, tavoutid, start, count, tancp_var);
00686   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00687     
00688 
00689   /* Close the output netCDF file */
00690   istat = ncclose(ncoutid);
00691   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00692     
00693 
00694 
00695   /* Open NetCDF file for writing, overwrite and truncate existing file if any */
00696   istat = nc_create(data->learning->filename_save_weight, NC_CLOBBER, &ncoutid);
00697   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00698 
00699   /* Set global attributes */
00700   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "processor", strlen(data->info->processor), data->info->processor);
00701   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "institution", strlen(data->info->institution), data->info->institution);
00702   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_email", strlen(data->info->creator_email), data->info->creator_email);
00703   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_url", strlen(data->info->creator_url), data->info->creator_url);
00704   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_name", strlen(data->info->creator_name), data->info->creator_name);
00705   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "contact_email", strlen(data->info->contact_email), data->info->contact_email);
00706   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "contact_name", strlen(data->info->contact_name), data->info->contact_name);
00707   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "other_contact_email", strlen(data->info->other_contact_email),
00708                           data->info->other_contact_email);
00709   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "other_contact_name", strlen(data->info->other_contact_name),
00710                           data->info->other_contact_name);
00711 
00712   fillvalue = -9999.9;
00713 
00714   /* Set dimensions */
00715   istat = nc_def_dim(ncoutid, "season", data->conf->nseasons, &sdimoutid);
00716   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00717   istat = nc_def_dim(ncoutid, data->conf->eofname, data->learning->rea_neof, &eofdimoutid);
00718   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00719 
00720   for (s=0; s<data->conf->nseasons; s++) {
00721 
00722     /* Define weight dimensions and variables */
00723     (void) sprintf(nomvar, "%s_%d", data->conf->clustname, s+1);
00724     istat = nc_def_dim(ncoutid, nomvar, data->conf->season[s].nclusters, &(weightdimoutid[s]));
00725     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00726 
00727     vardimids[0] = weightdimoutid[s];
00728     vardimids[1] = eofdimoutid;
00729     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_weight, s+1);
00730     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 2, vardimids, &(weightoutid[s]));
00731     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00732 
00733     istat = nc_put_att_double(ncoutid, weightoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00734     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00735     istat = sprintf(tmpstr, "%s %s_%d", data->conf->eofname, data->conf->clustname, s+1);
00736     istat = nc_put_att_text(ncoutid, weightoutid[s], "coordinates", strlen(tmpstr), tmpstr);
00737     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00738     istat = sprintf(tmpstr, "none");
00739     istat = nc_put_att_text(ncoutid, weightoutid[s], "units", strlen(tmpstr), tmpstr);
00740     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00741   }
00742 
00743   /* End definition mode */
00744   istat = nc_enddef(ncoutid);
00745   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00746 
00747   for (s=0; s<data->conf->nseasons; s++) {
00748     /* Write weights */
00749     start[0] = 0;
00750     start[1] = 0;
00751     start[2] = 0;
00752     count[0] = (size_t) data->conf->season[s].nclusters;
00753     count[1] = (size_t) data->learning->rea_neof;
00754     count[2] = 0;
00755     istat = nc_put_vara_double(ncoutid, weightoutid[s], start, count, data->learning->data[s].weight);
00756     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00757   }  
00758 
00759   /* Close the output netCDF file */
00760   istat = ncclose(ncoutid);
00761   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00762 
00763 
00764   /* Open NetCDF file for writing, overwrite and truncate existing file if any */
00765   istat = nc_create(data->learning->filename_save_clust_learn, NC_CLOBBER, &ncoutid);
00766   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00767 
00768   /* Set global attributes */
00769   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "processor", strlen(data->info->processor), data->info->processor);
00770   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "institution", strlen(data->info->institution), data->info->institution);
00771   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_email", strlen(data->info->creator_email), data->info->creator_email);
00772   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_url", strlen(data->info->creator_url), data->info->creator_url);
00773   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "creator_name", strlen(data->info->creator_name), data->info->creator_name);
00774   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "contact_email", strlen(data->info->contact_email), data->info->contact_email);
00775   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "contact_name", strlen(data->info->contact_name), data->info->contact_name);
00776   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "other_contact_email", strlen(data->info->other_contact_email),
00777                           data->info->other_contact_email);
00778   istat = nc_put_att_text(ncoutid, NC_GLOBAL, "other_contact_name", strlen(data->info->other_contact_name),
00779                           data->info->other_contact_name);
00780 
00781   fillvalue = -9999.9;
00782 
00783   /* Set dimensions */
00784   istat = nc_def_dim(ncoutid, "season", data->conf->nseasons, &sdimoutid);
00785   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00786   istat = nc_def_dim(ncoutid, data->conf->eofname, data->learning->rea_neof, &eofdimoutid);
00787   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00788 
00789   for (s=0; s<data->conf->nseasons; s++) {
00790 
00791     /* Define time dimensions and variables */
00792     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_time, s+1);
00793     istat = nc_def_dim(ncoutid, nomvar, data->learning->data[s].ntime, &(timedimoutid[s]));
00794     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00795 
00796     vardimids[0] = timedimoutid[s];
00797     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &(timeoutid[s]));
00798     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00799 
00800     istat = sprintf(tmpstr, "gregorian");
00801     istat = nc_put_att_text(ncoutid, timeoutid[s], "calendar", strlen(tmpstr), tmpstr);
00802     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00803     istat = sprintf(tmpstr, "%s", data->conf->time_units);
00804     istat = nc_put_att_text(ncoutid, timeoutid[s], "units", strlen(tmpstr), tmpstr);
00805     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00806     istat = sprintf(tmpstr, "time in %s", data->conf->time_units);
00807     istat = nc_put_att_text(ncoutid, timeoutid[s], "long_name", strlen(tmpstr), tmpstr);
00808     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00809 
00810     /* Define clust_learn variables */
00811     (void) sprintf(nomvar, "%s_%d", data->learning->nomvar_class_clusters, s+1);
00812     vardimids[0] = timedimoutid[s];
00813     istat = nc_def_var(ncoutid, nomvar, NC_DOUBLE, 1, vardimids, &(clustoutid[s]));
00814     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00815 
00816     istat = nc_put_att_double(ncoutid, clustoutid[s], "missing_value", NC_DOUBLE, 1, &fillvalue);
00817     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00818     istat = sprintf(tmpstr, "%s_%d", data->learning->nomvar_time, s+1);
00819     istat = nc_put_att_text(ncoutid, clustoutid[s], "coordinates", strlen(tmpstr), tmpstr);
00820     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00821     istat = sprintf(tmpstr, "none");
00822     istat = nc_put_att_text(ncoutid, clustoutid[s], "units", strlen(tmpstr), tmpstr);
00823     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00824   }
00825 
00826   /* End definition mode */
00827   istat = nc_enddef(ncoutid);
00828   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00829 
00830   for (s=0; s<data->conf->nseasons; s++) {
00831 
00832     timeval = (double *) realloc(timeval, data->learning->data[s].ntime * sizeof(double));
00833     if (timeval == NULL) alloc_error(__FILE__, __LINE__);
00834     
00835     /* Compute time variable using actual units */
00836     for (t=0; t<data->learning->data[s].ntime; t++)
00837       istat = utInvCalendar2(data->learning->data[s].time_s->year[t], data->learning->data[s].time_s->month[t],
00838                              data->learning->data[s].time_s->day[t], data->learning->data[s].time_s->hour[t],
00839                              data->learning->data[s].time_s->minutes[t], data->learning->data[s].time_s->seconds[t],
00840                              dataunits, &(timeval[t]));
00841 
00842     /* Write clust_learn */
00843     start[0] = 0;
00844     start[1] = 0;
00845     start[2] = 0;
00846     count[0] = (size_t) data->learning->data[s].ntime;
00847     count[1] = 0;
00848     count[2] = 0;
00849     istat = nc_put_vara_int(ncoutid, clustoutid[s], start, count, data->learning->data[s].class_clusters);
00850     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00851 
00852     /* Write time */
00853     start[0] = 0;
00854     start[1] = 0;
00855     start[2] = 0;
00856     count[0] = (size_t) data->learning->data[s].ntime;
00857     count[1] = 0;
00858     count[2] = 0;
00859     istat = nc_put_vara_double(ncoutid, timeoutid[s], start, count, data->learning->data[s].time);
00860     if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00861   }  
00862 
00863   /* Close the output netCDF file */
00864   istat = ncclose(ncoutid);
00865   if (istat != NC_NOERR) handle_netcdf_error(istat, __FILE__, __LINE__);
00866 
00867   (void) ut_free(dataunits);
00868   (void) ut_free_system(unitSystem);  
00869 
00870   (void) free(nomvar);
00871   (void) free(tancp_mean);
00872   (void) free(tancp_var);
00873   (void) free(timeval);
00874 
00875   (void) free(timedimoutid);
00876   (void) free(clustdimoutid);
00877   (void) free(weightdimoutid);
00878 
00879   (void) free(timeoutid);
00880   (void) free(cstoutid);
00881   (void) free(regoutid);
00882   (void) free(distoutid);
00883   (void) free(rrdoutid);
00884   (void) free(rrooutid);
00885   (void) free(taoutid);
00886   (void) free(tadoutid);
00887   (void) free(clustoutid);
00888   (void) free(weightoutid);
00889   (void) free(rsqoutid);
00890   (void) free(erroutid);
00891   (void) free(acoroutid);
00892   (void) free(vifoutid);
00893 
00894   (void) free(tmpstr);
00895 
00896   return 0;
00897 }

Generated on 12 May 2016 for DSCLIM by  doxygen 1.6.1