source: branches/ms-style/zoo-project/zoo-kernel/service_internal_ms.c @ 830

Last change on this file since 830 was 830, checked in by cresson, 7 years ago

ENH: Enhance raster styling in mapserver support

  • Property svn:eol-style set to native
  • Property svn:mime-type set to text/x-csrc
File size: 44.1 KB
Line 
1/*
2 * Author : Gérald FENOY
3 *
4 *  Copyright 2010-2011 Fondazione Edmund Mach. All rights reserved.
5 *
6 * Permission is hereby granted, free of charge, to any person obtaining a copy
7 * of this software and associated documentation files (the "Software"), to deal
8 * in the Software without restriction, including without limitation the rights
9 * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
10 * copies of the Software, and to permit persons to whom the Software is
11 * furnished to do so, subject to the following conditions:
12 *
13 * The above copyright notice and this permission notice shall be included in
14 * all copies or substantial portions of the Software.
15 *
16 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
17 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
18 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
19 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
20 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
21 * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
22 * THE SOFTWARE.
23 */
24
25/**
26 * Cross platform definition of layerObj->class
27 */
28#ifndef WIN32
29#define CLASS class
30#else
31#define CLASS _class
32#endif
33#include "service_internal_ms.h"
34#include "server_internal.h"
35#include "response_print.h"
36
37static double PRECISION = 0.001;
38static int MAX_NUMBER_STRING_SIZE = 32;
39
40/**
41 * Double to ASCII
42 */
43char * dtoa(char *s, double n) {
44
45  // handle special cases
46  if (isnan(n)) {
47    strcpy(s, "nan");
48  } else if (isinf(n)) {
49    strcpy(s, "inf");
50  } else if (n == 0.0) {
51    strcpy(s, "0");
52  } else {
53    int digit, m, m1;
54    char *c = s;
55    int neg = (n < 0);
56    if (neg)
57      n = -n;
58    // calculate magnitude
59    m = log10(n);
60    int useExp = (m >= 14 || (neg && m >= 9) || m <= -9);
61    if (neg)
62      *(c++) = '-';
63    // set up for scientific notation
64    if (useExp) {
65      if (m < 0)
66        m -= 1.0;
67      n = n / pow(10.0, m);
68      m1 = m;
69      m = 0;
70    }
71    if (m < 1.0) {
72      m = 0;
73    }
74    // convert the number
75    while (n > PRECISION || m >= 0) {
76      double weight = pow(10.0, m);
77      if (weight > 0 && !isinf(weight)) {
78        digit = floor(n / weight);
79        n -= (digit * weight);
80        *(c++) = '0' + digit;
81      }
82      if (m == 0 && n > 0)
83        *(c++) = '.';
84      m--;
85    }
86    if (useExp) {
87      // convert the exponent
88      int i, j;
89      *(c++) = 'e';
90      if (m1 > 0) {
91        *(c++) = '+';
92      } else {
93        *(c++) = '-';
94        m1 = -m1;
95      }
96      m = 0;
97      while (m1 > 0) {
98        *(c++) = '0' + m1 % 10;
99        m1 /= 10;
100        m++;
101      }
102      c -= m;
103      for (i = 0, j = m-1; i<j; i++, j--) {
104        // swap without temporary
105        c[i] ^= c[j];
106        c[j] ^= c[i];
107        c[i] ^= c[j];
108      }
109      c += m;
110    }
111    *(c) = '\0';
112  }
113  return s;
114}
115
116/**
117 * List of allowed raster styles
118 */
119enum MS_RASTER_STYLES{
120  LINEAR_STRETCHING,
121  COLOR_PALETTE
122};
123enum LINEAR_STRETCHING_TYPE{
124  MINMAX,
125  MEANSTD
126};
127
128/**
129 * Get a list of configuration keys having a corresponding mandatory ows_*.
130 * Map composed by a main.cfg maps name as key and the corresponding
131 * MapServer Mafile Metadata name to use
132 * see doc from here :
133 *  - http://mapserver.org/ogc/wms_server.html
134 *  - http://mapserver.org/ogc/wfs_server.html
135 *  - http://mapserver.org/ogc/wcs_server.html
136 *
137 * @return a new map containing a table linking a name of a configuration key
138 * to a corresponding mandatory ows_* keyword (ie. "fees" => "ows_fees").
139 */
140map* getCorrespondance(){
141  map* res=createMap("encoding","ows_encoding");
142  addToMap(res,"abstract","ows_abstract");
143  addToMap(res,"title","ows_title");
144  addToMap(res,"keywords","ows_keywordlist");
145  addToMap(res,"fees","ows_fees");
146  addToMap(res,"accessConstraints","ows_accessconstraints");
147  addToMap(res,"providerName","ows_attribution_title");
148  addToMap(res,"providerSite","ows_service_onlineresource");
149  addToMap(res,"individualName","ows_contactperson");
150  addToMap(res,"positionName","ows_contactposition");
151  addToMap(res,"providerName","ows_contactorganization");
152  addToMap(res,"role","ows_role");
153  addToMap(res,"addressType","ows_addresstype");
154  addToMap(res,"addressCity","ows_city");
155  addToMap(res,"addressDeliveryPoint","ows_address");
156  addToMap(res,"addressPostalCode","ows_postcode");
157  addToMap(res,"addressAdministrativeArea","ows_stateorprovince");
158  addToMap(res,"addressCountry","ows_country");
159  addToMap(res,"phoneVoice","ows_contactvoicetelephone");
160  addToMap(res,"phoneFacsimile","ows_contactfacsimiletelephone");
161  addToMap(res,"addressElectronicMailAddress","ows_contactelectronicmailaddress");
162  // Missing Madatory Informations
163  addToMap(res,"hoursOfService","ows_hoursofservice");
164  addToMap(res,"contactInstructions","ows_contactinstructions");
165  return res;
166}
167
168/**
169 * Add width and height keys to an output maps containing the maximum width
170 * and height for displaying the full data extent.
171 * Restriction to an image having a size of 640x480 (width * height)
172 *
173 * @param output
174 * @param minx the lower left x coordinate
175 * @param miny the lower left y coordinate
176 * @param maxx the upper right x coordinate
177 * @param maxy the upper right y coordinate
178 */
179void setMapSize(maps* output,double minx,double miny,double maxx,double maxy){
180  double maxWidth=640;
181  double maxHeight=480;
182  double deltaX=maxx-minx;
183  double deltaY=maxy-miny;
184  double qWidth;
185  qWidth=maxWidth/deltaX;
186  double qHeight;
187  qHeight=maxHeight/deltaY;
188#ifdef DEBUGMS
189  fprintf(stderr,"deltaX : %.15f \ndeltaY : %.15f\n",deltaX,deltaY);
190  fprintf(stderr,"qWidth : %.15f \nqHeight : %.15f\n",qWidth,qHeight);
191#endif
192
193  double width=deltaX*qWidth;
194  double height=height=deltaY*qWidth;
195  if(deltaX<deltaY){
196    width=deltaX*qHeight;
197    height=deltaY*qHeight;
198  }
199  if(height<0)
200    height=-height;
201  if(width<0)
202    width=-width;
203  char sWidth[1024];
204  char sHeight[1024];
205  sprintf(sWidth,"%.3f",width);
206  sprintf(sHeight,"%.3f",height);
207#ifdef DEBUGMS
208  fprintf(stderr,"sWidth : %.15f \nsHeight : %.15f\n",sWidth,sHeight);
209#endif
210  if(output!=NULL){
211    addToMap(output->content,"width",sWidth);
212    addToMap(output->content,"height",sHeight);
213  }
214}
215
216/**
217 * Add a Reference key to an output containing the WMFS/WFS/WCS request for
218 * accessing service result
219 *
220 * @param m the conf maps containing the main.cfg settings
221 * @param tmpI the specific output maps to add the Reference key
222 */
223void setReferenceUrl(maps* m,maps* tmpI){
224  outputMapfile(m,tmpI);
225  map *msUrl=getMapFromMaps(m,"main","mapserverAddress");
226  if(msUrl==NULL){
227    errorException (m, _("Unable to find any mapserverAddress defined in the main.cfg file"),
228                    "InternalError", NULL);
229    exit(-1);
230  }
231  map *msOgcVersion=getMapFromMaps(m,"main","msOgcVersion");
232  map *dataPath=getMapFromMaps(m,"main","dataPath");
233  map *sid=getMapFromMaps(m,"lenv","usid");
234  map* format=getMap(tmpI->content,"mimeType");
235  map* rformat=getMap(tmpI->content,"requestedMimeType");
236  map* width=getMap(tmpI->content,"width");
237  map* height=getMap(tmpI->content,"height");
238  map* protoMap=getMap(tmpI->content,"msOgc");
239  map* versionMap=getMap(tmpI->content,"msOgcVersion");
240  char options[3][5][25]={
241    {"WMS","1.3.0","GetMap","layers=%s","wms_extent"},
242    {"WFS","1.1.0","GetFeature","typename=%s","wcs_extent"},
243    {"WCS","1.1.0","GetCoverage","coverage=%s","wcs_extent"}
244  };
245  int proto=0;
246  if(rformat==NULL){
247    rformat=getMap(tmpI->content,"mimeType");
248  }
249  if(strncasecmp(rformat->value,"text/xml",8)==0)
250    proto=1;
251  if(strncasecmp(rformat->value,"image/tiff",10)==0)
252    proto=2;
253  if(protoMap!=NULL){
254    if(strncasecmp(protoMap->value,"WMS",3)==0)
255      proto=0;
256    else{
257      if(strncasecmp(protoMap->value,"WFS",3)==0)
258        proto=1;
259      else 
260        proto=2;
261    }
262  }
263  char *protoVersion=options[proto][1];
264  if(proto==1){
265    if(msOgcVersion!=NULL)
266      protoVersion=msOgcVersion->value;
267    if(versionMap!=NULL)
268      protoVersion=versionMap->value;
269  }
270
271  map* extent=getMap(tmpI->content,options[proto][4]);
272  map* crs=getMap(tmpI->content,"crs");
273  int hasCRS=1;
274  if(crs==NULL){
275    crs=getMapFromMaps(m,"main","crs");
276    if(crs==NULL){
277      crs=createMap("crs","epsg:4326");
278      hasCRS=0;
279    }
280  }
281  char layers[128];
282  sprintf(layers,options[proto][3],tmpI->name);
283
284  char* webService_url=(char*)malloc((strlen(msUrl->value)+strlen(format->value)+strlen(tmpI->name)+strlen(width->value)+strlen(height->value)+strlen(extent->value)+256)*sizeof(char));
285
286  if(proto>0){
287    sprintf(webService_url,
288            "%s?map=%s/%s_%s.map&request=%s&service=%s&version=%s&%s&format=%s&bbox=%s&crs=%s",
289            msUrl->value,
290            dataPath->value,
291            tmpI->name,
292            sid->value,
293            options[proto][2],
294            options[proto][0],
295            protoVersion,
296            layers,
297            rformat->value,
298            extent->value,
299            crs->value
300            );
301  }
302  else{
303    sprintf(webService_url,
304            "%s?map=%s/%s_%s.map&request=%s&service=%s&version=%s&%s&width=%s&height=%s&format=%s&bbox=%s&crs=%s",
305            msUrl->value,
306            dataPath->value,
307            tmpI->name,
308            sid->value,
309            options[proto][2],
310            options[proto][0],
311            protoVersion,
312            layers,
313            width->value,
314            height->value,
315            rformat->value,
316            extent->value,
317            crs->value
318            );
319  }
320  if(hasCRS==0){
321    freeMap(&crs);
322    free(crs);
323  }
324  addToMap(tmpI->content,"Reference",webService_url);
325  free(webService_url);
326}
327
328/**
329 * Set projection for a layer in a MAPFILE using Authority Code and Name if
330 * available or fallback to proj4 definition if available or fallback to
331 * default EPSG:4326
332 *
333 * @param output the output maps
334 * @param m the opened mapObj
335 * @param myLayer the layerObj
336 * @param pszProjection a char* containing the SRS definition in WKT format
337 */
338void setSrsInformations(maps* output,mapObj* m,layerObj* myLayer,
339                        char* pszProjection){
340  OGRSpatialReferenceH  hSRS;
341  map* msSrs=NULL;
342  hSRS = OSRNewSpatialReference(NULL);
343  if( pszProjection!=NULL && strlen(pszProjection)>1){
344    if(OSRImportFromWkt( hSRS, &pszProjection ) == CE_None ){
345      char *proj4Str=NULL;
346      if(OSRGetAuthorityName(hSRS,NULL)!=NULL && 
347         OSRGetAuthorityCode(hSRS,NULL)!=NULL){
348        char tmpSrs[20];
349        sprintf(tmpSrs,"%s:%s",
350                OSRGetAuthorityName(hSRS,NULL),OSRGetAuthorityCode(hSRS,NULL));
351        msLoadProjectionStringEPSG(&m->projection,tmpSrs);
352        msLoadProjectionStringEPSG(&myLayer->projection,tmpSrs);
353       
354        char tmpSrss[256];
355        sprintf(tmpSrss,"EPSG:4326 EPSG:900913 EPSG:3857 %s",tmpSrs);
356       
357        msInsertHashTable(&(m->web.metadata), "ows_srs", tmpSrss);
358        msInsertHashTable(&(myLayer->metadata), "ows_srs", tmpSrss);
359       
360#ifdef DEBUGMS
361        fprintf(stderr,"isGeo %b\n\n",OSRIsGeographic(hSRS)==TRUE);
362#endif
363        if(output!=NULL){
364          if(OSRIsGeographic(hSRS)==TRUE)
365            addToMap(output->content,"crs_isGeographic","true");
366          else
367            addToMap(output->content,"crs_isGeographic","false");
368          addToMap(output->content,"crs",tmpSrs);
369        }
370      }
371      else{
372        OSRExportToProj4(hSRS,&proj4Str);
373        if(proj4Str!=NULL && strlen(proj4Str)>0){
374#ifdef DEBUGMS
375          fprintf(stderr,"PROJ (%s)\n",proj4Str);
376#endif
377          msLoadProjectionString(&(m->projection),proj4Str);
378          msLoadProjectionString(&(myLayer->projection),proj4Str);
379          if(output!=NULL){ 
380            if(OSRIsGeographic(hSRS)==TRUE)
381              addToMap(output->content,"crs_isGeographic","true");
382            else
383              addToMap(output->content,"crs_isGeographic","false");
384          }
385          free(proj4Str);
386        }
387        else{
388          msLoadProjectionStringEPSG(&m->projection,"EPSG:4326");
389          msLoadProjectionStringEPSG(&myLayer->projection,"EPSG:4326");
390          if(output!=NULL){
391            addToMap(output->content,"crs_isGeographic","true");
392          }
393        }
394        if(output!=NULL){
395          addToMap(output->content,"crs","EPSG:4326");
396          addToMap(output->content,"real_extent","true");
397        }
398        msInsertHashTable(&(m->web.metadata),"ows_srs", "EPSG:4326 EPSG:900913 EPSG:3857");
399        msInsertHashTable(&(myLayer->metadata),"ows_srs","EPSG:4326 EPSG:900913 EPSG:3857");
400      }
401    }
402  }
403  else{
404    if(output!=NULL){
405      msSrs=getMap(output->content,"msSrs");
406    }
407    if(msSrs!=NULL){
408      msLoadProjectionStringEPSG(&m->projection,msSrs->value);
409      msLoadProjectionStringEPSG(&myLayer->projection,msSrs->value);
410      char tmpSrs[128];
411      sprintf(tmpSrs,"%s EPSG:4326 EPSG:900913 EPSG:3857",msSrs->value);
412      msInsertHashTable(&(m->web.metadata),"ows_srs",tmpSrs);
413      msInsertHashTable(&(myLayer->metadata),"ows_srs",tmpSrs);
414    }else{
415      msLoadProjectionStringEPSG(&m->projection,"EPSG:4326");
416      msLoadProjectionStringEPSG(&myLayer->projection,"EPSG:4326");
417      msInsertHashTable(&(m->web.metadata),"ows_srs","EPSG:4326 EPSG:900913 EPSG:3857");
418      msInsertHashTable(&(myLayer->metadata),"ows_srs","EPSG:4326 EPSG:900913 EPSG:3857");
419    }
420    if(output!=NULL){
421      addToMap(output->content,"crs",msSrs->value);
422      addToMap(output->content,"crs_isGeographic","true");
423    }
424  }
425
426  OSRDestroySpatialReference( hSRS );
427}
428
429/**
430 * Set the MAPFILE extent, the the ows_extent for the layer, add wms_extent and
431 * wfs_extent to the output maps and call setMapSize.
432 *
433 * @param output the specific output
434 * @param m the mapObj
435 * @param myLayer the layerObj
436 * @param minX the lower left x coordinate
437 * @param minY the lower left y coordinate
438 * @param maxX the upper right x coordinate
439 * @param maxY the upper right y coordinate
440 * @see setMapSize
441 */
442void setMsExtent(maps* output,mapObj* m,layerObj* myLayer,
443                 double minX,double minY,double maxX,double maxY){
444  msMapSetExtent(m,minX,minY,maxX,maxY);
445#ifdef DEBUGMS
446  fprintf(stderr,"Extent %.15f %.15f %.15f %.15f\n",minX,minY,maxX,maxY);
447#endif
448  char tmpExtent[1024];
449  sprintf(tmpExtent,"%.15f %.15f %.15f %.15f",minX,minY,maxX,maxY);
450#ifdef DEBUGMS
451  fprintf(stderr,"Extent %s\n",tmpExtent);
452#endif
453  msInsertHashTable(&(myLayer->metadata), "ows_extent", tmpExtent);
454 
455  if(output!=NULL){
456    map* test=getMap(output->content,"real_extent");
457    if(test!=NULL){
458      pointObj min, max;
459      projectionObj tempSrs;
460      min.x = m->extent.minx;
461      min.y = m->extent.miny;
462      max.x = m->extent.maxx;
463      max.y = m->extent.maxy;
464      char tmpSrsStr[1024];
465      msInitProjection(&tempSrs);
466      msLoadProjectionStringEPSG(&tempSrs,"EPSG:4326");
467
468      msProjectPoint(&(m->projection),&tempSrs,&min);
469      msProjectPoint(&m->projection,&tempSrs,&max);
470     
471      sprintf(tmpExtent,"%.3f,%.3f,%.3f,%.3f",min.y,min.x,max.y,max.x);
472      map* isGeo=getMap(output->content,"crs_isGeographic");
473#ifdef DEBUGMS
474      fprintf(stderr,"isGeo = %s\n",isGeo->value);
475#endif
476      if(isGeo!=NULL && strcasecmp("true",isGeo->value)==0)
477        sprintf(tmpExtent,"%f,%f,%f,%f", minY,minX, maxY, maxX);
478      addToMap(output->content,"wms_extent",tmpExtent);
479      sprintf(tmpSrsStr,"%.3f,%.3f,%.3f,%.3f",min.x,min.y,max.x,max.y);
480      addToMap(output->content,"wcs_extent",tmpExtent);
481    }else{
482      sprintf(tmpExtent,"%f,%f,%f,%f",minX, minY, maxX, maxY);
483      map* isGeo=getMap(output->content,"crs_isGeographic");
484      if(isGeo!=NULL){
485#ifdef DEBUGMS
486        fprintf(stderr,"isGeo = %s\n",isGeo->value);
487#endif
488        if(isGeo!=NULL && strcasecmp("true",isGeo->value)==0)
489          sprintf(tmpExtent,"%f,%f,%f,%f", minY,minX, maxY, maxX);
490      }
491      addToMap(output->content,"wms_extent",tmpExtent); 
492      sprintf(tmpExtent,"%.3f,%.3f,%.3f,%.3f",minX,minY,maxX,maxY);
493      addToMap(output->content,"wcs_extent",tmpExtent);
494    }
495  }
496
497  setMapSize(output,minX,minY,maxX,maxY);
498}
499
500/**
501 * Try to open a vector output and define the corresponding layer in the MAPFILE
502 *
503 * @param conf the conf maps containing the main.cfg settings
504 * @param output the specific output maps
505 * @param m the mapObj
506 */
507int tryOgr(maps* conf,maps* output,mapObj* m){
508
509  map* tmpMap=getMap(output->content,"storage");
510  char *pszDataSource=tmpMap->value;
511
512  /**
513   * Try to open the DataSource using OGR
514   */
515  OGRRegisterAll();
516  /**
517   * Try to load the file as ZIP
518   */
519
520  OGRDataSourceH poDS1 = NULL;
521  OGRSFDriverH *poDriver1 = NULL;
522  char *dsName=(char*)malloc((8+strlen(pszDataSource)+1)*sizeof(char));
523  char *odsName=zStrdup(pszDataSource);
524  char *sdsName=zStrdup(pszDataSource);
525  char *demo=".data";
526  sdsName[strlen(sdsName)-(strlen(demo)-1)]='d';
527  sdsName[strlen(sdsName)-(strlen(demo)-2)]='i';
528  sdsName[strlen(sdsName)-(strlen(demo)-3)]='r';
529  sdsName[strlen(sdsName)-(strlen(demo)-4)]=0;
530
531  odsName[strlen(odsName)-(strlen(demo)-1)]='z';
532  odsName[strlen(odsName)-(strlen(demo)-2)]='i';
533  odsName[strlen(odsName)-(strlen(demo)-3)]='p';
534  odsName[strlen(odsName)-(strlen(demo)-4)]=0;
535  sprintf(dsName,"/vsizip/%s",odsName);
536
537#ifdef DEBUGMS
538  fprintf(stderr,"Try loading %s, %s, %s\n",dsName,odsName,dsName);
539#endif
540
541  FILE* file = fopen(pszDataSource, "rb");
542  FILE* fileZ = fopen(odsName, "wb");
543  free(odsName);
544  fseek(file, 0, SEEK_END);
545  unsigned long fileLen=ftell(file);
546  fseek(file, 0, SEEK_SET);
547  char *buffer=(char *)malloc(fileLen+1);
548  fread(buffer, fileLen, 1, file);
549  fwrite(buffer,fileLen, 1, fileZ);
550  fclose(file);
551  fclose(fileZ);
552  free(buffer);
553#ifdef DEBUGMS
554  fprintf(stderr,"Try loading %s",dsName);
555#endif
556  poDS1 = OGROpen( dsName, FALSE, poDriver1 );
557  if( poDS1 == NULL ){
558    fprintf(stderr,"Unable to access the DataSource as ZIP File\n");
559    setMapInMaps(conf,"lenv","message","Unable to open datasource in read only mode");
560    OGR_DS_Destroy(poDS1);
561  }else{
562#ifdef DEBUGMS
563    fprintf(stderr,"The DataSource is a  ZIP File\n");
564#endif
565    char** demo=VSIReadDir(dsName);
566    int i=0;
567    zMkdir(sdsName
568#ifndef WIN32
569                ,S_IRWXU | S_IRGRP | S_IXGRP | S_IROTH | S_IXOTH
570#endif
571                );
572    while(demo[i]!=NULL){
573#ifdef DEBUGMS
574      fprintf(stderr,"ZIP File content : %s\n",demo[i]);
575#endif
576      char *tmpDs=(char*)malloc((strlen(dsName)+strlen(demo[i])+2)*sizeof(char));
577      sprintf(tmpDs,"%s/%s",dsName,demo[i]);
578      fprintf(stderr,"read : %s\n",tmpDs);
579     
580      VSILFILE* vsif=VSIFOpenL(tmpDs,"rb");
581#ifdef DEBUGMS
582      fprintf(stderr,"open : %s\n",tmpDs);
583#endif
584      VSIFSeekL(vsif,0,SEEK_END);
585      vsi_l_offset size=VSIFTellL(vsif);
586#ifdef DEBUGMS
587      fprintf(stderr,"size : %d\n",size);
588#endif
589      VSIFSeekL(vsif,0,SEEK_SET);
590      char *vsifcontent=(char*) malloc(((int)size+1)*sizeof(char));
591      VSIFReadL(vsifcontent,1,(size_t)size,vsif);
592      char *fpath=(char*) malloc((strlen(sdsName)+strlen(demo[1])+2)*sizeof(char));
593      sprintf(fpath,"%s/%s",sdsName,demo[i]);
594      int f=zOpen(fpath,O_WRONLY|O_CREAT,S_IRUSR|S_IWUSR|S_IRGRP|S_IWGRP|S_IROTH|S_IWOTH);
595      zWrite(f,vsifcontent,(int)size);
596      close(f);
597      chmod(fpath,S_IRWXU | S_IRGRP | S_IXGRP | S_IROTH);
598      char* tmpP=strstr(fpath,".shp");
599      if(tmpP==NULL)
600        tmpP=strstr(fpath,".SHP");
601      if(tmpP!=NULL){
602#ifdef DEBUGMS
603        fprintf(stderr,"*** DEBUG %s\n",strstr(tmpP,"."));
604#endif
605        if( strcmp(tmpP,".shp")==0 || strcmp(tmpP,".SHP")==0 ){
606          tmpMap=getMap(output->content,"storage");
607          free(tmpMap->value);
608          tmpMap->value=(char*) malloc((strlen(fpath)+1)*sizeof(char));
609          sprintf(tmpMap->value,"%s",fpath);
610          pszDataSource=tmpMap->value;
611#ifdef DEBUGMS
612          fprintf(stderr,"*** DEBUG %s\n",pszDataSource);
613#endif
614        }
615      }
616      VSIFCloseL(vsif);
617      i++;
618    }
619  }
620  free(sdsName);
621  free(dsName);
622
623  OGRDataSourceH poDS = NULL;
624  OGRSFDriverH *poDriver = NULL;
625  poDS = OGROpen( pszDataSource, FALSE, poDriver );
626  if( poDS == NULL ){
627#ifdef DEBUGMS
628    fprintf(stderr,"Unable to access the DataSource %s\n",pszDataSource);
629#endif
630    setMapInMaps(conf,"lenv","message","Unable to open datasource in read only mode");
631    OGR_DS_Destroy(poDS);
632    OGRCleanupAll();
633#ifdef DEBUGMS
634    fprintf(stderr,"Unable to access the DataSource, exit! \n"); 
635#endif
636    return -1;
637  }
638
639  int iLayer = 0;
640  for( iLayer=0; iLayer < OGR_DS_GetLayerCount(poDS); iLayer++ ){
641    OGRLayerH poLayer = OGR_DS_GetLayer(poDS,iLayer);
642
643    if( poLayer == NULL ){
644#ifdef DEBUGMS
645      fprintf(stderr,"Unable to access the DataSource Layer \n");
646#endif
647      setMapInMaps(conf,"lenv","message","Unable to open datasource in read only mode");
648      return -1;
649    }
650
651    /**
652     * Add a new layer set name, data
653     */
654    if(msGrowMapLayers(m)==NULL){
655      return -1;
656    }
657    if(initLayer((m->layers[m->numlayers]), m) == -1){
658      return -1;
659    }
660
661    layerObj* myLayer=m->layers[m->numlayers];
662#ifdef DEBUGMS
663    dumpMaps(output);
664#endif
665    myLayer->name = zStrdup(output->name);
666    myLayer->tileitem=NULL;
667    myLayer->data = zStrdup(OGR_L_GetName(poLayer));
668    myLayer->connection = zStrdup(pszDataSource);
669    myLayer->index = m->numlayers;
670    myLayer->dump = MS_TRUE;
671    myLayer->status = MS_ON;
672    msConnectLayer(myLayer,MS_OGR,pszDataSource);
673
674    /**
675     * Detect the Geometry Type or use Polygon
676     */
677    if(OGR_L_GetGeomType(poLayer) != wkbUnknown){
678      switch(OGR_L_GetGeomType(poLayer)){
679      case wkbPoint:
680      case wkbMultiPoint:
681      case wkbPoint25D:
682      case wkbMultiPoint25D:
683#ifdef DEBUGMS
684        fprintf(stderr,"POINT DataSource Layer \n");
685#endif
686        myLayer->type = MS_LAYER_POINT;
687        break;
688      case wkbLineString :
689      case wkbMultiLineString :
690      case wkbLineString25D:
691      case wkbMultiLineString25D:
692#ifdef DEBUGMS
693        fprintf(stderr,"LINE DataSource Layer \n");
694#endif
695        myLayer->type = MS_LAYER_LINE;
696        break;
697      case wkbPolygon:
698      case wkbMultiPolygon:
699      case wkbPolygon25D:
700      case wkbMultiPolygon25D:
701#ifdef DEBUGMS
702        fprintf(stderr,"POLYGON DataSource Layer \n");
703#endif
704        myLayer->type = MS_LAYER_POLYGON;
705        break;
706      default:
707        myLayer->type = MS_LAYER_POLYGON;
708        break;
709      }
710    }else
711      myLayer->type = MS_LAYER_POLYGON;
712
713    /**
714     * Detect spatial reference or use WGS84
715     */
716    OGRSpatialReferenceH srs=OGR_L_GetSpatialRef(poLayer);
717    if(srs!=NULL){
718      char *wkt=NULL;
719      OSRExportToWkt(srs,&wkt);
720      setSrsInformations(output,m,myLayer,wkt);
721      free(wkt);
722    }
723    else{
724      addToMap(output->content,"crs","EPSG:4326");
725      addToMap(output->content,"crs_isGeographic","true");
726      msLoadProjectionStringEPSG(&m->projection,"EPSG:4326");
727      msInsertHashTable(&(m->web.metadata), "ows_srs", "EPSG:4326 EPSG:900913 EPSG:3857");
728      msInsertHashTable(&(myLayer->metadata), "ows_srs", "EPSG:4326 EPSG:900913 EPSG:3857");
729    }
730
731    OGREnvelope oExt;
732    if (OGR_L_GetExtent(poLayer,&oExt, TRUE) == OGRERR_NONE){
733      setMsExtent(output,m,myLayer,oExt.MinX, oExt.MinY, oExt.MaxX, oExt.MaxY);
734    }
735 
736    /**
737     * Detect the FID column or use the first attribute field as FID
738     */
739    char *fid=(char*)OGR_L_GetFIDColumn(poLayer);
740    if(strlen(fid)==0){
741      OGRFeatureDefnH def=OGR_L_GetLayerDefn(poLayer);
742      int fIndex=0;
743      for(fIndex=0;fIndex<OGR_FD_GetFieldCount(def);fIndex++){
744        OGRFieldDefnH fdef=OGR_FD_GetFieldDefn(def,fIndex);
745        fid=(char*)OGR_Fld_GetNameRef(fdef);
746        break;
747      }
748    }
749    msInsertHashTable(&(myLayer->metadata), "gml_featureid", fid);
750    msInsertHashTable(&(myLayer->metadata), "gml_include_items", "all");
751    msInsertHashTable(&(myLayer->metadata), "ows_name", output->name);
752    map* tmpMap=getMap(output->content,"title");
753    if(tmpMap!=NULL)
754      msInsertHashTable(&(myLayer->metadata), "ows_title", tmpMap->value);
755    else
756      msInsertHashTable(&(myLayer->metadata), "ows_title", "Default Title");
757   
758    if(msGrowLayerClasses(myLayer) == NULL)
759      return -1;
760    if(initClass((myLayer->CLASS[myLayer->numclasses])) == -1)
761      return -1;
762    if(msGrowClassStyles(myLayer->CLASS[myLayer->numclasses]) == NULL)
763      return -1;
764    if(initStyle(myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]) == -1)
765      return -1;
766    /**
767     * Apply msStyle else fallback to the default style
768     */
769    tmpMap=getMap(output->content,"msStyle");
770    if(tmpMap!=NULL)
771      msUpdateStyleFromString(myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles],tmpMap->value,0);
772    else{
773      /**
774       * Set style
775       */
776      myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->color.red=125;
777      myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->color.green=125;
778      myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->color.blue=255;
779      myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->outlinecolor.red=80;
780      myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->outlinecolor.green=80;
781      myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->outlinecolor.blue=80;
782     
783      /**
784       * Set specific style depending on type
785       */
786      if(myLayer->type == MS_LAYER_POLYGON)
787        myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->width=3;
788      if(myLayer->type == MS_LAYER_LINE){
789        myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->width=3;
790        myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->outlinewidth=1.5;
791      }
792      if(myLayer->type == MS_LAYER_POINT){
793        myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->symbol=1;
794        myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->size=15;
795      }
796     
797    }
798    myLayer->CLASS[myLayer->numclasses]->numstyles++;
799    myLayer->numclasses++;
800   
801    m->layerorder[m->numlayers] = m->numlayers;
802    m->numlayers++;
803
804  }
805
806  OGR_DS_Destroy(poDS);
807  OGRCleanupAll();
808
809  return 1;
810}
811
812/**
813 * Try to open a raster output and define the corresponding layer in the MAPFILE
814 *
815 * @param conf the conf maps containing the main.cfg settings
816 * @param output the specific output maps
817 * @param m the mapObj
818 */
819int tryGdal(maps* conf,maps* output,mapObj* m){
820
821  /*
822   * Detect the raster style
823   */
824
825  /* msRasterResample (NEAREST/AVERAGE/BILINEAR) */
826  const char * msRasterResamplingType               = "msRasterResample";
827  /* msRasterStyle (linearStretching/colorPalette) */
828  const char * msRasterStylePropertyName            = "msRasterStyle";
829  const char * msRasterStyleLinearStretchingPropertyValue       = "linearStretching";
830  const char * msRasterStyleColorPalettePropertyValue           = "colorPalette";
831  /* msRasterStyleOptions (minMax/meanstd) */
832  const char * msRasterStyleOptionsPropertyName     = "msRasterStyleOptions";
833  const char * msRasterStyleLinearStretchingMinMaxPropertyName  = "minMax";
834  const char * msRasterStyleLinearStretchingMeanStdPropertyName = "meanStd";
835
836  // Default raster style
837  int defaultStyleType = LINEAR_STRETCHING;
838  int defaultLinearstretchingType = MEANSTD;
839
840  // Check if there is a defined raster style type
841  int styleType = defaultStyleType;
842  int linearStretchingType = defaultLinearstretchingType;
843  map* msRasterStyle=getMap(output->content, msRasterStylePropertyName);
844  char * msRasterStyleOptionsContent = "";
845  if(msRasterStyle!=NULL)
846    {
847    // Check if there is options attached
848    map* msRasterStyleOptions=getMap(output->content, msRasterStyleOptionsPropertyName);
849    if (msRasterStyleOptions!=NULL)
850      {
851      msRasterStyleOptionsContent = msRasterStyleOptions->value;
852      }
853
854    // Detect the raster style
855#ifdef DEBUGMS
856    fprintf(stderr,"Detect the raster style %s\n", msRasterStyle->value);
857#endif
858    if (strncasecmp(msRasterStyle->value, msRasterStyleLinearStretchingPropertyValue,
859        strlen(msRasterStyleLinearStretchingPropertyValue))==0)
860      {
861#ifdef DEBUGMS
862      fprintf(stderr,"The raster style is linear stretching\n");
863#endif
864      styleType = LINEAR_STRETCHING;
865      if (strlen(msRasterStyleOptionsContent)>0)
866        {
867        if (strncasecmp(msRasterStyleOptionsContent, msRasterStyleLinearStretchingMinMaxPropertyName,
868            strlen(msRasterStyleLinearStretchingMinMaxPropertyName))==0)
869          {
870          linearStretchingType = MINMAX;
871#ifdef DEBUGMS
872          fprintf(stderr,"The raster style linear stretching is minmax\n");
873#endif
874          }
875        else if (strncasecmp(msRasterStyleOptionsContent, msRasterStyleLinearStretchingMeanStdPropertyName,
876            strlen(msRasterStyleLinearStretchingMeanStdPropertyName))==0)
877          {
878          linearStretchingType = MEANSTD;
879#ifdef DEBUGMS
880          fprintf(stderr,"The raster style linear stretching is meanstd\n");
881#endif
882          }
883        else
884          {
885          fprintf(stderr,"Unknown raster style linear stretching method: %s\n", msRasterStyleOptionsContent);
886          }
887        }
888      else
889        {
890        fprintf(stderr,"Using default linear stretching type.");
891        }
892      }
893    else if (strncasecmp(msRasterStyle->value, msRasterStyleColorPalettePropertyValue,
894        strlen(msRasterStyleColorPalettePropertyValue))==0)
895      {
896#ifdef DEBUGMS
897      fprintf(stderr,"The raster style is color palette\n");
898#endif
899      styleType = COLOR_PALETTE;
900      }
901    else
902      {
903      fprintf(stderr,"Unknown raster style: %s. Using default.", msRasterStyle->value);
904      }
905    }
906
907#ifdef DEBUGMS
908  fprintf(stderr,"RasterStyle=%i Options=%s\n",styleType,msRasterStyleOptionsContent);
909#endif
910
911  /*
912   * Get storage
913   */
914  map* tmpMap=getMap(output->content,"storage");
915  char *pszFilename=tmpMap->value;
916  GDALDatasetH hDataset;
917  GDALRasterBandH hBand;
918  double adfGeoTransform[6];
919  int i, iBand;
920
921  /**
922   * Try to open the DataSource using GDAL
923   */
924  GDALAllRegister();
925  hDataset = GDALOpen( pszFilename, GA_ReadOnly );
926  if( hDataset == NULL ){
927#ifdef DEBUGMS
928    fprintf(stderr,"Unable to access the DataSource %s \n",pszFilename);
929#endif
930    setMapInMaps(conf,"lenv","message","gdalinfo failed - unable to open");
931    GDALDestroyDriverManager();
932    return -1;
933  }
934#ifdef DEBUGMS
935  fprintf(stderr,"Accessing the DataSource %s %d\n",pszFilename,__LINE__);
936#endif
937
938  /**
939   * Add a new layer set name, data
940   */
941  if(msGrowMapLayers(m)==NULL){
942    return -1;
943  }
944  if(initLayer((m->layers[m->numlayers]), m) == -1){
945    return -1;
946  }
947  m->layers[m->numlayers]->index=m->numlayers;
948
949  layerObj* myLayer=m->layers[m->numlayers];
950  myLayer->name = zStrdup(output->name);
951  myLayer->tileitem=NULL;
952  myLayer->data = zStrdup(pszFilename);
953  myLayer->index = m->numlayers;
954  myLayer->dump = MS_TRUE;
955  myLayer->status = MS_ON;
956  myLayer->type = MS_LAYER_RASTER;
957
958  char *title=output->name;
959  tmpMap=getMap(output->content,"title");
960  if(tmpMap!=NULL)
961    title=tmpMap->value;
962  char *abstract=output->name;
963  tmpMap=getMap(output->content,"abstract");
964  if(tmpMap!=NULL)
965    abstract=tmpMap->value;
966
967  msInsertHashTable(&(myLayer->metadata), "ows_label", title);
968  msInsertHashTable(&(myLayer->metadata), "ows_title", title);
969  msInsertHashTable(&(myLayer->metadata), "ows_abstract", abstract);
970  msInsertHashTable(&(myLayer->metadata), "ows_rangeset_name", output->name);
971  msInsertHashTable(&(myLayer->metadata), "ows_rangeset_label", title);
972
973  /**
974   * Set Map Size to the raster size
975   */
976  m->width=GDALGetRasterXSize( hDataset );
977  m->height=GDALGetRasterYSize( hDataset );
978
979  /**
980   * Set projection using Authority Code and Name if available or fallback to
981   * proj4 definition if available or fallback to default EPSG:4326
982   */
983  const char *tRef=GDALGetProjectionRef( hDataset );
984  if( tRef != NULL && strlen(tRef)>0 ){
985    char *pszProjection;
986    pszProjection = (char *) GDALGetProjectionRef( hDataset );
987#ifdef DEBUGMS
988    fprintf(stderr,"Accessing the DataSource %s\n",GDALGetProjectionRef( hDataset ));
989#endif
990    setSrsInformations(output,m,myLayer,pszProjection);
991  }else{
992    fprintf(stderr,"NO SRS FOUND ! %s\n",GDALGetProjectionRef( hDataset ));   
993  }
994
995  /**
996   * Set extent
997   */
998  if( GDALGetGeoTransform( hDataset, adfGeoTransform ) == CE_None ){
999    if( adfGeoTransform[2] == 0.0 && adfGeoTransform[4] == 0.0 ){
1000
1001      double minX = adfGeoTransform[0]
1002                                    + adfGeoTransform[2] * GDALGetRasterYSize(hDataset);
1003      double minY = adfGeoTransform[3]
1004                                    + adfGeoTransform[5] * GDALGetRasterYSize(hDataset);
1005
1006      double maxX = adfGeoTransform[0]
1007                                    + adfGeoTransform[1] * GDALGetRasterXSize(hDataset);
1008      double maxY = adfGeoTransform[3]
1009                                    + adfGeoTransform[4] * GDALGetRasterXSize(hDataset);
1010
1011      setMsExtent(output,m,myLayer,minX,minY,maxX,maxY);
1012
1013    }
1014  }
1015
1016  /**
1017   * Extract information about available bands to set the bandcount and the
1018   * processing directive
1019   */
1020  char nBands[3];
1021  int nBandsI=GDALGetRasterCount( hDataset );
1022  sprintf(nBands,"%d",GDALGetRasterCount( hDataset ));
1023  msInsertHashTable(&(myLayer->metadata), "ows_bandcount", nBands);
1024  if(nBandsI>=3)
1025    msLayerAddProcessing(myLayer,"BANDS=1,2,3");
1026  else if(nBandsI>=2)
1027    msLayerAddProcessing(myLayer,"BANDS=1,2");
1028  else
1029    msLayerAddProcessing(myLayer,"BANDS=1");
1030
1031  /**
1032   * Name available Bands
1033   */
1034  char lBands[7];
1035  char *nameBands=NULL;
1036  for( iBand = 0; iBand < nBandsI; iBand++ ){
1037    sprintf(lBands,"Band%d",iBand+1);
1038    if(nameBands==NULL){
1039      nameBands=(char*)malloc((strlen(lBands)+1)*sizeof(char));
1040      sprintf(nameBands,"%s",lBands);
1041    }else{
1042      if(iBand<4){
1043        char *tmpS=zStrdup(nameBands);
1044        nameBands=(char*)realloc(nameBands,(strlen(nameBands)+strlen(lBands)+1)*sizeof(char));
1045        sprintf(nameBands,"%s %s",tmpS,lBands);
1046        free(tmpS);
1047      }
1048    }
1049  }
1050  msInsertHashTable(&(myLayer->metadata), "ows_bandnames", nameBands);
1051
1052  /**
1053   * Loops over metadata information to setup specific information
1054   */
1055  for( iBand = 0; iBand < nBandsI; iBand++ )
1056    {
1057
1058    // Compute statistics of the current band
1059    hBand = GDALGetRasterBand( hDataset, iBand+1 );
1060    CPLErrorReset();
1061    double min = 0.0;
1062    double max = 0.0;
1063    double mean = 0.0;
1064    double std = 0.0;
1065    GDALComputeRasterStatistics  (hBand, 1, &min, &max, &mean, &std, NULL, NULL);
1066
1067    char bandIdentifier[21];
1068    sprintf(bandIdentifier,"Band%d",iBand+1);
1069    if (CPLGetLastErrorType() == CE_None)
1070      {
1071      char bandInterval[100];
1072      sprintf(bandInterval,"%.3f %.3f",min,max);
1073      char bandIntervalIdentifier[21];
1074      sprintf(bandIntervalIdentifier,"%s_interval",bandIdentifier);
1075      msInsertHashTable(&(myLayer->metadata), bandIntervalIdentifier, bandInterval);
1076
1077      // Apply the raster style
1078      if(styleType == LINEAR_STRETCHING)
1079        {
1080
1081        char msProcessingInstruction[1024];
1082        double low = mean - 2*std;
1083        double hi = mean + 2*std;
1084
1085        char s1[MAX_NUMBER_STRING_SIZE];
1086        char s2[MAX_NUMBER_STRING_SIZE];
1087
1088        int bn = iBand+1;
1089        if(linearStretchingType==MINMAX)
1090          {
1091          sprintf(msProcessingInstruction, "SCALE_%d=%s,%s",
1092              bn,
1093              dtoa(s1,min),
1094              dtoa(s2,max));
1095          }
1096        else if (linearStretchingType==MEANSTD)
1097          {
1098          sprintf(msProcessingInstruction, "SCALE_%d=%s,%s",
1099              bn,
1100              dtoa(s1,low),
1101              dtoa(s2,hi));
1102          }
1103
1104        msLayerAddProcessing(myLayer,msProcessingInstruction);
1105
1106        } // styleType is LINEAR_STRETCHING
1107      else if( styleType == COLOR_PALETTE )
1108        {
1109        /**
1110         * Classify one band raster pixel value using regular interval
1111         * TODO: parse msRasterStyleOptionsContent to retrieve the colormap
1112         */
1113        int discreteColorMap[10][3]={
1114            {102,153,204},
1115            {51,102,153},
1116            {102,102,204},
1117            {51,204,0},
1118            {153,255,102},
1119            {204,255,102},
1120            {102,204,153},
1121            {255,69,64},
1122            {255,192,115},
1123            {255,201,115}
1124        };
1125
1126        if(nBandsI==1)
1127          {
1128          double delta = max - min;
1129          double interval = delta / 10;
1130          double cstep = min;
1131          for(i=0; i<10; i++)
1132            {
1133            /**
1134             * Create a new class
1135             */
1136            if(msGrowLayerClasses(myLayer) == NULL)
1137              return -1;
1138            if(initClass((myLayer->CLASS[myLayer->numclasses])) == -1)
1139              return -1;
1140            if(msGrowClassStyles(myLayer->CLASS[myLayer->numclasses]) == NULL)
1141              return -1;
1142            if(initStyle(myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]) == -1)
1143              return -1;
1144
1145            /**
1146             * Set class name
1147             */
1148            char className[7];
1149            sprintf(className,"class%d",i);
1150            myLayer->CLASS[myLayer->numclasses]->name=zStrdup(className);
1151
1152            /**
1153             * Set expression
1154             */
1155            char expression[1024];
1156            if(i+1<10)
1157              sprintf(expression,"([pixel]>=%.3f AND [pixel]<%.3f)",cstep,cstep+interval);
1158            else
1159              sprintf(expression,"([pixel]>=%.3f AND [pixel]<=%.3f)",cstep,cstep+interval);
1160            msLoadExpressionString(&myLayer->CLASS[myLayer->numclasses]->expression,expression);
1161
1162            /**
1163             * Set color
1164             */
1165            myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->color.red=discreteColorMap[i][0];
1166            myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->color.green=discreteColorMap[i][1];
1167            myLayer->CLASS[myLayer->numclasses]->styles[myLayer->CLASS[myLayer->numclasses]->numstyles]->color.blue=discreteColorMap[i][2];
1168            cstep+=interval;
1169            myLayer->CLASS[myLayer->numclasses]->numstyles++;
1170            myLayer->numclasses++;
1171
1172            }
1173
1174          char tmpMm[100];
1175          sprintf(tmpMm,"%.3f %.3f",min,max);
1176
1177          }
1178        } // styleType is COLOR_PALETTE
1179
1180      } // If no error with GDAL functions
1181    else
1182      {
1183      fprintf(stderr,"Unable to compute raster statistics!\n");
1184      }
1185
1186    /*
1187     * Set offsite
1188     */
1189    int offsiteR = 0;
1190    int offsiteG = 0;
1191    int offsiteB = 0;
1192    int hasNoData = 0;
1193    double noDataValue = GDALGetRasterNoDataValue(hBand, &hasNoData);
1194    if (hasNoData)
1195      {
1196      offsiteR = (int) noDataValue;
1197      offsiteG = (int) noDataValue;
1198      offsiteB = (int) noDataValue;
1199      }
1200    myLayer->offsite.red    = offsiteR;
1201    myLayer->offsite.green  = offsiteG;
1202    myLayer->offsite.blue   = offsiteB;
1203
1204    /*
1205     * Insert units
1206     */
1207    if( strlen(GDALGetRasterUnitType(hBand)) > 0 )
1208      {
1209      char tmpU[21];
1210      sprintf(tmpU,"%s_band_uom",bandIdentifier);
1211      msInsertHashTable(&(myLayer->metadata), tmpU, GDALGetRasterUnitType(hBand));
1212      }
1213
1214    } // next band
1215
1216  /*
1217   * Check if there is resample option
1218   */
1219  char msResampleOptionInstruction[1024];
1220  char * msRasterResampleOptionContent = "BILINEAR";
1221  map* msRasterResamplingOption=getMap(output->content, msRasterResamplingType);
1222  if (msRasterResamplingOption!=NULL)
1223    {
1224    msRasterResampleOptionContent = msRasterResamplingOption->value;
1225    }
1226  sprintf(msResampleOptionInstruction, "RESAMPLE=%s",msRasterResampleOptionContent);
1227  msLayerAddProcessing(myLayer, msResampleOptionInstruction);
1228
1229  m->layerorder[m->numlayers] = m->numlayers;
1230  m->numlayers++;
1231  GDALClose( hDataset );
1232  GDALDestroyDriverManager();
1233  CPLCleanupTLS();
1234  return 1;
1235}
1236
1237/**
1238 * Create a MapFile for WMS, WFS or WCS Service output
1239 *
1240 * @param conf the conf maps containing the main.cfg settings
1241 * @param outputs a specific output maps
1242 */
1243void outputMapfile(maps* conf,maps* outputs){
1244
1245  /**
1246   * First store the value on disk
1247   */
1248  map* mime=getMap(outputs->content,"mimeType");
1249  char *ext="data";
1250  if(mime!=NULL)
1251    if(strncasecmp(mime->value,"application/json",16)==0)
1252      ext="json";
1253
1254  map* tmpMap=getMapFromMaps(conf,"main","dataPath");
1255  map* sidMap=getMapFromMaps(conf,"lenv","usid");
1256  char *pszDataSource=(char*)malloc((strlen(tmpMap->value)+strlen(sidMap->value)+strlen(outputs->name)+17)*sizeof(char));
1257  sprintf(pszDataSource,"%s/ZOO_DATA_%s_%s.%s",tmpMap->value,outputs->name,sidMap->value,ext); 
1258  int f=zOpen(pszDataSource,O_WRONLY|O_CREAT,S_IRUSR|S_IWUSR|S_IRGRP|S_IWGRP|S_IROTH|S_IWOTH);
1259  map *gfile=getMap(outputs->content,"generated_file");
1260  if(gfile!=NULL){
1261    readGeneratedFile(conf,outputs->content,gfile->value);         
1262  }
1263  map* sizeMap=getMap(outputs->content,"size");
1264  map* vData=getMap(outputs->content,"value");
1265  if(sizeMap!=NULL){
1266    zWrite(f,vData->value,atoi(sizeMap->value)*sizeof(char));
1267  }
1268  else{
1269    zWrite(f,vData->value,(strlen(vData->value)+1)*sizeof(char));
1270  }
1271  close(f);
1272  addToMap(outputs->content,"storage",pszDataSource);
1273  free(pszDataSource);
1274
1275  /*
1276   * Create an empty map, set name, default size and extent
1277   */
1278  mapObj *myMap=msNewMapObj();
1279  free(myMap->name);
1280  myMap->name=zStrdup("ZOO-Project_WXS_Server");
1281  msMapSetSize(myMap,2048,2048);
1282  msMapSetExtent(myMap,-1,-1,1,1);
1283 
1284  /*
1285   * Set imagepath and imageurl using tmpPath and tmpUrl from main.cfg
1286   */
1287  map *tmp1=getMapFromMaps(conf,"main","tmpPath");
1288  myMap->web.imagepath=zStrdup(tmp1->value);
1289  tmp1=getMapFromMaps(conf,"main","tmpUrl");
1290  myMap->web.imageurl=zStrdup(tmp1->value);
1291 
1292  /*
1293   * Define supported output formats
1294   */
1295  outputFormatObj *o1=msCreateDefaultOutputFormat(NULL,"AGG/PNG","png");
1296  o1->imagemode=MS_IMAGEMODE_RGBA;
1297  o1->transparent=MS_TRUE;
1298  o1->inmapfile=MS_TRUE;
1299  msAppendOutputFormat(myMap,msCloneOutputFormat(o1));
1300  msFreeOutputFormat(o1);
1301
1302#ifdef USE_KML
1303  outputFormatObj *o2=msCreateDefaultOutputFormat(NULL,"KML","kml");
1304  o2->inmapfile=MS_TRUE; 
1305  msAppendOutputFormat(myMap,msCloneOutputFormat(o2));
1306  msFreeOutputFormat(o2);
1307#endif
1308
1309  outputFormatObj *o3=msCreateDefaultOutputFormat(NULL,"GDAL/GTiff","tiff");
1310  if(!o3)
1311    fprintf(stderr,"Unable to initialize GDAL driver !\n");
1312  else{
1313    o3->imagemode=MS_IMAGEMODE_BYTE;
1314    o3->inmapfile=MS_TRUE; 
1315    msAppendOutputFormat(myMap,msCloneOutputFormat(o3));
1316    msFreeOutputFormat(o3);
1317  }
1318
1319  outputFormatObj *o4=msCreateDefaultOutputFormat(NULL,"GDAL/AAIGRID","grd");
1320  if(!o4)
1321    fprintf(stderr,"Unable to initialize GDAL driver !\n");
1322  else{
1323    o4->imagemode=MS_IMAGEMODE_INT16;
1324    o4->inmapfile=MS_TRUE; 
1325    msAppendOutputFormat(myMap,msCloneOutputFormat(o4));
1326    msFreeOutputFormat(o4);
1327  }
1328
1329#ifdef USE_CAIRO
1330  outputFormatObj *o5=msCreateDefaultOutputFormat(NULL,"CAIRO/PNG","cairopng");
1331  if(!o5)
1332    fprintf(stderr,"Unable to initialize CAIRO driver !\n");
1333  else{
1334    o5->imagemode=MS_IMAGEMODE_RGBA;
1335    o5->transparent=MS_TRUE;
1336    o5->inmapfile=MS_TRUE;
1337    msAppendOutputFormat(myMap,msCloneOutputFormat(o5));
1338    msFreeOutputFormat(o5);
1339  }
1340#endif
1341
1342  /*
1343   * Set default projection to EPSG:4326
1344   */
1345  msLoadProjectionStringEPSG(&myMap->projection,"EPSG:4326");
1346  myMap->transparent=1;
1347
1348  /**
1349   * Set metadata extracted from main.cfg file maps
1350   */
1351  maps* cursor=conf;
1352  map* correspondance=getCorrespondance();
1353  while(cursor!=NULL){
1354    map* _cursor=cursor->content;
1355    map* vMap;
1356    while(_cursor!=NULL){
1357      if((vMap=getMap(correspondance,_cursor->name))!=NULL){
1358        if (msInsertHashTable(&(myMap->web.metadata), vMap->value, _cursor->value) == NULL){
1359#ifdef DEBUGMS
1360          fprintf(stderr,"Unable to add metadata");
1361#endif
1362          return;
1363        }
1364      }
1365      _cursor=_cursor->next;
1366    }
1367    cursor=cursor->next;
1368  }
1369  freeMap(&correspondance);
1370  free(correspondance);
1371
1372  /*
1373   * Set mapserver PROJ_LIB/GDAL_DATA or any other config parameter from
1374   * the main.cfg [mapserver] section if any
1375   */
1376  maps *tmp3=getMaps(conf,"mapserver");
1377  if(tmp3!=NULL){
1378    map* tmp4=tmp3->content;
1379    while(tmp4!=NULL){
1380      msSetConfigOption(myMap,tmp4->name,tmp4->value);
1381      tmp4=tmp4->next;
1382    }
1383  }
1384
1385  /**
1386   * Set a ows_rootlayer_title, 
1387   */
1388  if (msInsertHashTable(&(myMap->web.metadata), "ows_rootlayer_name", "ZOO_Project_Layer") == NULL){
1389#ifdef DEBUGMS
1390    fprintf(stderr,"Unable to add metadata");
1391#endif
1392    return;
1393  }
1394  if (msInsertHashTable(&(myMap->web.metadata), "ows_rootlayer_title", "ZOO_Project_Layer") == NULL){
1395#ifdef DEBUGMS
1396    fprintf(stderr,"Unable to add metadata");
1397#endif
1398    return;
1399  }
1400
1401  /**
1402   * Enable all the WXS requests using ows_enable_request
1403   * see http://mapserver.org/trunk/development/rfc/ms-rfc-67.html
1404   */
1405  if (msInsertHashTable(&(myMap->web.metadata), "ows_enable_request", "*") == NULL){
1406#ifdef DEBUGMS
1407    fprintf(stderr,"Unable to add metadata");
1408#endif
1409    return;
1410  }
1411  msInsertHashTable(&(myMap->web.metadata), "ows_srs", "EPSG:4326");
1412
1413  if(tryOgr(conf,outputs,myMap)<0)
1414    if(tryGdal(conf,outputs,myMap)<0)
1415      return ;
1416
1417  tmp1=getMapFromMaps(conf,"main","dataPath");
1418  char *tmpPath=(char*)malloc((13+strlen(tmp1->value))*sizeof(char));
1419  sprintf(tmpPath,"%s/symbols.sym",tmp1->value);
1420  msInitSymbolSet(&myMap->symbolset);
1421  myMap->symbolset.filename=zStrdup(tmpPath);
1422  free(tmpPath);
1423
1424  map* sid=getMapFromMaps(conf,"lenv","usid");
1425  char *mapPath=
1426    (char*)malloc((7+strlen(sid->value)+strlen(outputs->name)+strlen(tmp1->value))*sizeof(char));
1427  sprintf(mapPath,"%s/%s_%s.map",tmp1->value,outputs->name,sid->value);
1428  msSaveMap(myMap,mapPath);
1429  free(mapPath);
1430  msGDALCleanup();
1431  msFreeMap(myMap);
1432}
1433
Note: See TracBrowser for help on using the repository browser.

Search

ZOO Sponsors

http://www.zoo-project.org/trac/chrome/site/img/geolabs-logo.pnghttp://www.zoo-project.org/trac/chrome/site/img/neogeo-logo.png http://www.zoo-project.org/trac/chrome/site/img/apptech-logo.png http://www.zoo-project.org/trac/chrome/site/img/3liz-logo.png http://www.zoo-project.org/trac/chrome/site/img/gateway-logo.png

Become a sponsor !

Knowledge partners

http://www.zoo-project.org/trac/chrome/site/img/ocu-logo.png http://www.zoo-project.org/trac/chrome/site/img/gucas-logo.png http://www.zoo-project.org/trac/chrome/site/img/polimi-logo.png http://www.zoo-project.org/trac/chrome/site/img/fem-logo.png http://www.zoo-project.org/trac/chrome/site/img/supsi-logo.png http://www.zoo-project.org/trac/chrome/site/img/cumtb-logo.png

Become a knowledge partner

Related links

http://zoo-project.org/img/ogclogo.png http://zoo-project.org/img/osgeologo.png