{'fraction_non_alphanumeric': 0.047490771226179784, 'fraction_numerical': 0.03912667531344574, 'mean_word_length': 4.119084099421178, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 13, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
NGC 7538 : Multiwavelength Study of Stellar Cluster Regions associated with IRS 1-3 and IRS 9 sources

We present deep and high-resolution (FWHM ~ 0.4 arcsec) near-infrared (NIR) imaging observations of the NGC 7538 IRS 1-3 region (in JHK bands), and IRS 9 region (in HK bands) using the 8.2m Subaru telescope. The NIR analysis is complemented with GMRT low-frequency observations at 325, 610, and 1280 MHz, molecular line observations of H13CO+ (J=1-0), and archival Chandra X-ray observations. Using the 'J-H/H-K' diagram, 144 Class II and 24 Class I young stellar object (YSO) candidates are identified in the IRS 1-3 region. Further analysis using 'K/H-K' diagram yields 145 and 96 red sources in the IRS 1-3 and IRS 9 regions, respectively. A total of 27 sources are found to have X-ray counterparts. The YSO mass function (MF), constructed using a theoretical mass-luminosity relation, shows peaks at substellar (~0.08-0.18 Msolar) and intermediate (~1-1.78 Msolar) mass ranges for the IRS 1-3 region. The MF can be fitted by a power law in the low mass regime with a slope of Gamma ~ 0.54-0.75, which is much shallower than the Salpeter value of 1.35. An upper limit of 10.2 is obtained for the star to brown dwarf ratio in the IRS 1-3 region. GMRT maps show a compact HII region associated with the IRS 1-3 sources, whose spectral index of 0.87+-0.11 suggests optical thickness. This compact region is resolved into three separate peaks in higher resolution 1280 MHz map, and the 'East' sub-peak coincides with the IRS 2 source. H13CO+ (J=1-0) emission reveals peaks in both IRS 1-3 and IRS 9 regions, none of which are coincident with visible nebular emission, suggesting the presence of dense cloud nearby. The virial masses are approximately of the order of 1000 Msolar and 500 Msolar for the clumps in IRS 1-3 and IRS 9 regions, respectively.

detail using various techniques such as high-frequency radio observations, line emission observations, maser observations, outflows, and so on. However, most of these studies as well as many recent ones concentrate on the individual features of luminous IRS and candidate high-mass sources. Even then, the bulk deals with the sources around IRS 1-3 stellar cluster region, and the IRS 9 region analyses have been sparse. Recently, an NIR study of the overall NGC 7538 region by Ojha et al. (2004a), reddening and cluster related studies using NIR by Balog et al. (2004), a spectroscopic study of the luminous sources by Puga et al. (2010), and a multiwavelength study of clusters (using statistical techniques) with a focus on high-mass stars by Chavarría et al. (2014) have been carried out. Other large scale studies have been carried out at far-infrared Herschel bands and at submillimetre wavelengths to identify cold clumps and filamentary structures (Sandell & Sievers 2004;Fallscheer et al. 2013).
According to Ojha et al. (2004a) and McCaughrean et al. (1991), there appear to be three distinct sub-regions which can be separated based on the stellar population and morphology, namely the IRS 1-3, IRS 4-8, and IRS 9 subregions. Although the previous works identified a rich cluster membership in a wide field, encompassing all the major starforming sites in this complex, none of the optical/infrared surveys were deep enough to reach the substellar regime. Hence the works available in the literature were not able to have a detailed study of these stellar cluster sub-regions, which need to be analysed separately due to their veritable differences. Recent advanced instruments on 8-10 m class telescopes now make it feasible to conduct imaging and spectroscopic studies of low-mass populations in distant highmass star-forming regions, where the star-forming environment may be different from their low-mass counterparts. Moreover, in the case of stars in a distant cluster, high resolution imaging is also needed to recognize individual stars. Since most of the available studies of low-mass populations are mainly based on nearby star-forming regions and the sample of substellar sources is small to draw definitive conclusions, we considered it worthwhile to focus on the stellar clusters of the NGC 7538 region observed using the 8.2 m Subaru telescope -the deepest and highest-resolution NIR data till date. The NIR data has been used to study the luminosity function and initial mass function (IMF) of the region. In young massive star-forming regions, a gamut of components are present which could affect further evolution, and thus we have complemented our deep NIR observations with : previously unexamined low-frequency bremsstrahlung to ascertain the ionizing gas characteristics, molecular line H 13 CO + (J=1 -0) observations for the dense gas morphology, and Chandra X-ray observations for stellar population analysis.
In this paper, therefore, the aim is to continue with our multiwavelength study of star-forming regions (Ojha et al. 2004a(Ojha et al. , 2011Samal et al. 2007Samal et al. , 2010, as well as our investigations to detect and characterize the young brown dwarfs (BDs) in distant massive star-forming regions (cf. Ojha et al. 2004bOjha et al. , 2009. In Section 2, we present details of the observations and data reduction procedures. The YSO selection procedure is dealt with in Section 3. The spatial distribution, morphology, and physical characteristics of the regions are discussed in Section 4. In Section 5, we elaborate upon the luminosity and mass functions obtained for differ-ent clusters. Discussion and final conclusions are presented in Sections 6 and 7, respectively.

Near Infrared Photometry
Deep NIR imaging observations of the NGC 7538 IRS 1-3 region (centered on α2000 = 23 h 13 m 43 s , δ2000 = +61 o 28 22 ) in J (λ =1.25 µm), H (λ =1.64 µm), and K (λ =2.21 µm) bands, and the NGC 7538 IRS 9 region (centered on α2000 = 23 h 13 m 58 s , δ2000 = +61 o 27 26 ) in H and K bands were obtained on 2005 August 19, using the Cooled Infrared Spectrograph and Camera for OHS (CISCO) mounted at the Cassegrain focus of the 8.2 m Subaru telescope. These observations were done in the service mode of the telescope. CISCO is equipped with a 1024×1024 Rockwell HgCdTe HAWAII array. A plate scale of 0.105 pixel −1 at the f/12 focus of the telescope provides a field-of-view (FoV) of ∼ 1.8 × 1.8 (Motohara et al. 2002). Observations were carried out in a 3×3 dithering pattern with ∼ 10 offsets -with a varying number of images (three upto ten) being obtained at each dithered position. Individual image exposure times were 40 s, 20 s, and 10 s for the J, H, and K bands, respectively, finally yielding a total integration time of 12, 12, and 13.5 minutes in these respective bands. All the observations were done under excellent photometric sky conditions. The average seeing size was measured to be 0.4 full-width-athalf-maxima (FWHM) in all three filters, and the air mass variation was between 1.34 and 1.44. Off-target (located ∼ 34 east of the target position) sky-frame observations, identical in area to the target FoV, were taken just after the target observations using a similar procedure.
Data reduction was done using the Image Reduction and Analysis Facility (iraf) software package. The sky flats were generated by median-combining individual dithered sky frames for respective filters. These median-combined sky-flats were applied for both flat-fielding and sky subtraction. Identification of the point sources and their photometry was performed using the daofind and daophot packages of iraf. Because of source confusion and nebulosity within the region, photometry was performed using the point spread function (PSF) algorithm allstar in the daophot package (Stetson 1987). An aperture radius of 4 pixels (∼ 0.4 ) was used for the final photometry, with appropriate aperture corrections applied for the respective bands. Since no standard star was observed during the observations, the photometric calibration was carried out using sources (about 10-12, cross-matched within 0.4 of the Subaru catalogue) from the Two Micron All Sky Survey (2MASS) 1 , selected on the basis of their 'phqual' and 'ccflg' flag values, as well as visual examination, to avoid artifacts and contaminations. Finally, the photometric calibration rms was ∼ 0.1 mag. The Subaru/CISCO system magnitudes are assumed to be in California Institute of Technology (CIT) system (Oasa et al. 2006), and hence the 2MASS magnitudes too were converted to CIT system (using Carpenter 2001) for calibration purpose. On comparison of our photometry with that of Ojha et al. (2004a), we find that the average dispersion (in the entire K magnitude range) was ∼ 0.05-0.10 mag for the JHK bands. Our higher spatial resolution permits better source separation and sky determination. Sources which were bright and saturated in Subaru images, but had good quality ('phqual=A or B' for J, H, and Ks bands each) 2MASS photometry had their magnitudes taken from the 2MASS catalog after conversion to CIT system. Similar photometric procedure was also carried out for the off-target sky region. However, the western edge of the sky region was found to be slightly contaminated by a nearby nebula, and hence the sky frame was trimmed to remove the nebulous western portion before doing the photometry. The final sky frame size on which photometry (and completeness calculation below) was carried out is ∼ 1.17 ×1.8 . Since these observations were carried out in service mode, it was found that during the later observations, of IRS 9 region in K band, source profiles at the north-east part of the images were elongated, due to some likely optics problem. Since this part (north-east of IRS 9 region) of the image is not crowded or nebulous, aperture photometry was carried out to estimate the source magnitudes of the sources with elongated profiles. Even then, due to limitations, we use IRS 9 photomtery for qualitative assessments only in this work.
The completeness limits of the images were evaluated through artificial star experiments using addstar in iraf. Since the IRS 1-3 and IRS 9 regions have varying nebulosity, the images were divided into separate regions as shown in Figure 1, followed by completeness determination for each sub-region. A fixed number of stars were added in every 0.5 magnitude interval, followed by photometry to see how many of these added stars were recovered. This cycle was carried out repeatedly. We thus obtained the detection rate -which is just the ratio of the number of recovered artificial stars to the number of added stars -as a function of magnitude for each of the sub-regions in IRS 1-3 and IRS 9 regions, as well as the sky region. Table 1 summarizes the 90% completeness limits for all three bands in each of the sub-regions. The 10σ limiting magnitudes for our observations are estimated to be ∼ 22, 21 and 20 in the J, H, and K bands, respectively. As can be seen, the sky frame 90% completeness limits are equal to these 10 σ limiting magnitudes for the H and K bands, and slightly lower for the J band.

Radio Continuum Observations
Radio continuum observations were carried out using the Giant Metrewave Radio Telescope (GMRT) for the frequency bands 325 MHz (2004 July 03), 610 MHz (2004 September 18), and 1280 MHz (2004 January 25). The GMRT array, consisting of 30 antennae, is in an approximate Yshaped configuration. Each of these antennae has a diameter of 45 m, and thus a primary beam size of ∼ 81 , 43 , and 26.2 for 325, 610, and 1280 MHz, respectively 2 . There is a central region (∼ 1 km×1 km) which consists of randomly distributed 12 antennae, while the rest of the antennae are along three radial arms (6 along each arm) extending upto ∼ 14 km each. The minimum and maximum baselines are 100 m and 25 km, respectively. Further details about the GMRT array can be looked up in Swarup et al. (1991).
The total observation time ranged from ∼ 2.5-3.5 hours for the three frequency bands. Data reduction was carried out using the aips software. Successive rounds of flagging and calibration were carried out to improve the calibration. Flagging involved removal of bad data (including bad antennae, baselines, channels, time-ranges resulting from terrestrial radio frequency interference, etc), and was done using the 'vplot-uvflg' and 'tvflg' tasks. The respective flux and phase calibrators were used in the standard 'calib-getjy-clcal' procedure for amplitude and phase calibration. After satisfactory calibration, the source data (NGC 7538) was 'split' from the whole file (which contains flux and phase calibrator data in addition). Facet imaging was done using the aips task 'imagr' to generate the requisite images. To remove the ionospheric phase distortion effects, a few rounds of (phase) self-calibration -with decreasing 'solint' -were carried out using the task 'calib' . Table 2 gives the details of the GMRT observations and the parameters of the generated images. In addition, VLA archival image for 4860 MHz, for the observation date 2000 September 22 (Project ID BP0068), was also used in the analysis 3 .

Chandra X-ray Observations
Publicly available archival X-ray data for NGC 7538 region was retrieved from the Chandra site 4 (Obs. ID 5373). The X-ray observations of this H ii region were carried out using the Advanced CCD Imaging Spectrometer (ACIS; Garmire et al. 2003) onboard the Chandra X-ray Observatory (CXO; Weisskopf et al. 2002). For our purpose, we use only the imaging array of ACIS (ACIS-I). ACIS-I consists of a 2×2 CCD array of 1024×1024 pixels each, with a total FoV of ∼ 17 ×17 . The net exposure time was 30 ks.
Initial steps in the data reduction were carried out using the Chandra Interactive Analysis of Observations (ciao; Fruscione et al. 2006) tool version 4.5 and Chandra Calibration Database (caldb) version 4.5.8. The data was reprocessed using the 'chandra repro' tool to apply the latest calibration to it. Light curve of source-free background regions were constructed to verify that the data was not affected by solar flare activity. Subsequently, with the help of the wavelet based source detection tool -'wavdetect' (Freeman et al. 2002) -source detection was carried out at a threshold level of 2×10 −6 . Data image and exposure map with the default pixel scale of ∼ 0.5 were used for this purpose. Once we have a final list of positions of detected sources from 'wavdetect' , we used the IDL based ACIS Extract (ae) software package 5 , 2012) version 'March 6, 2013 to extract the relevant source properties. The algorithm detailed in the ae user's guide was followed. The source counts are extracted within ∼ 0.90 PSF fraction, with the PSF being separate for each source. Within ae, the X-ray spectra are compiled and fitted with an optically-thin thermal plasma model attenuated by an interstellar absorption using the 'xspec' fitting package. A total of 182 sources were obtained in the FoV, all of which were verified to have PB 0.01, to make sure that the detected sources are real and not a result of Poissonian fluctuation in the local background. PB gives the probability that the extracted counts in the total band are solely a result of background fluctuations. The intrinsic hard band luminosity for all the sources was found to range from ∼ 5.5×10 28 to 3.5×10 32 erg s −1 , while the intrinsic total luminosity was in the range ∼ 5.5×10 29 to 6×10 35 erg s −1 . The hard band luminosity function peak was at about 5×10 30 erg s −1 . For sources with spectral fitting results, the column density peaks at ∼ 10 22 cm −2 and the plasma temperature at ∼ 1 keV. The brightest ones show some metallic emission lines, supporting the thermal origin of the X-rays. These typical features are commonly seen among YSOs (Getman et al. 2006). We additionally used the non-parametric xphot program of Getman et al. (2010) to estimate the intrinsic fluxes and X-ray column densities.
Among the total number of sources obtained, there will also be extragalactic contaminants like Active Galactic Nuclei (AGN), and foreground stellar sources. Though we use a small subset of sources for our study (see Sections 3.1 and 3.2), and thus have not carried out a detailed contaminantion analysis here, we can obtain an estimate using the values for the Cepheus B region (Getman et al. 2006) which is close to NGC 7538 and had the same exposure time as well as FoV. Getman et al. (2006) find that the extragalactic contamination in Cepheus B is 5%, while stellar contamination is 4%. Assuming these numbers would imply a total of ∼ 16 contaminants (out of 182 total) in the NGC 7538 ACIS-I FoV. Table 3 gives the properties of 27 X-ray sources with NIR counterparts in the NIR FoV which have been selected in this paper for further analysis.

H 13 CO + (J=1 -0) Observations
The H 13 CO + (J=1-0) (formylium) molecular line (86.754 GHz) observations were carried out on 2004 May 02 with the Nobeyama 45 m radio telescope. At 87 GHz, the halfpower beam width and main beam efficiency, η, of the telescope were 18 and 0.51, respectively. We used the 25-BEam Array Receiver System (BEARS) (Sunada et al. 2000). To correct for the beam-to-beam gain variation, we calibrated the intensity scale of each beam using a 100 GHz SIS receiver (S100) with an SSB filter. Furthermore, we observed the same grid point in mapping with 9 different beams to smooth the beam-to-beam gain variation. At the back end, we used 25 arrays of 1024 channel Auto-Correlators (ACs), which have a 32 MHz band width and a 37.8 kHz resolution, corresponding to 0.13 km s −1 (Sorai et al. 2000). All the observations were carried out in position-switching mode. The standard chopper wheel method was used to convert the received intensity into the antenna temperature, T * A . Our mapping observations covered the same region of the NIR image. The mapping grid has 21 spacing, corresponding to half of the beam separation of the BEARS, i.e., nearly fullbeam sampling. During the observations, the system noise temperatures were in the range of 200 to 400 K, resulting in a noise level of 0.35 K in T * A . The telescope pointing was checked every 1.5 hours by observing the SiO maser source R Cas. The pointing accuracy was better than 3 .

Using NIR Colour-Colour Diagram
The photometric catalogs containing the J, H, and K band magnitudes for the IRS 1-3 region (see Section 2.1) and the sky field region were used to construct the NIR J-H vs H-K colour-colour diagrams (CC-Ds) shown in Figure 2. In this figure, the red curve denotes the dwarf locus from Bessell & Brett (1988) (converted to CIT system), blue solid line denotes the locus of Classical T Tauri Stars (CTTS) from Meyer et al. (1997), and the three parallel dashed green lines denote the reddening vectors drawn using the reddening laws of Cohen et al. (1981) (AJ /AV = 0.265, AH /AV = 0.155, and AK /AV = 0.090) for the CIT system. Three separate regions are marked on the images -'F' , 'T' , and 'P' -similar to Ojha et al. (2004a,b). The 'F' region mostly contains the field stars and Class III-type (Weak-line T Tauri or WTTS) sources, the 'T' region mostly contains the CTTS and Class II-type sources (Lada & Adams 1992), and the 'P' region mostly contains the Class I-type sources with circumstellar envelopes. Since the sources in the 'P' region can be contaminated by Herbig Ae-Be stars (Lada & Adams 1992), we conservatively only consider those sources in this region whose J-H colour is larger than that of the CTTS locus extended into this region. There may also be an overlap in the NIR colours of the upper end band of Herbig Ae-Be stars and in the lower end band of CTTS in the 'T' region (Hillenbrand et al. 1992).
The sky field CC-D shown in Figure 2(b) is used to examine the extent upto which IRS 1-3 CC-D is affected by field star contamination. As can clearly be seen in Figure  2(b), almost all field star contamination is present in the 'F' region. Therefore, while the sources falling in the 'T' and 'P' regions of Figure 2(a) are most probably YSOs, the Class IIItype sources from the 'F' region will be contaminated by field stars.
To further separate the most-likely YSOs from the 'F' region, since pre-main sequence (PMS) stars display much stronger X-ray emission than the field main sequence (MS) stars (Feigelson & Montmerle 1999;Montmerle 1996), we use the X-ray source identifications from Section 2.3. For this, we cross-matched our NIR catalog with the X-ray catalog within 0.5 radius. X-ray detected sources, however, in general, suffer from extragalactic contamination (mostly AGN) and foreground stellar contamination. But since we are using only those X-ray sources which have NIR counterparts, it is unlikely -similar to the work of Getman et al. (2006, for Cepheus B) and Kuhn et al. (2010, for W40) -that there will be any extragalactic contamination. Also since the resolution of Chandra images (∼ 0.5 ) is similar to our NIR imaging (∼ 0.4 , see Section 2.1), mismatching should not be a concern. From the matching radius and the source number density, the number of X-ray sources to have an NIR counterpart by chance was calculated to be one at most. Apropos foreground contamination, the X-ray sources lie beyond the low-density gap at about H − K ∼ 0.5 (again similar to Kuhn et al. 2010) in the NIR diagrams. Normally, a sudden low-density gap signifies the boundary of foreground sources in H − K colour space. Hence, these X-ray matched sources are unlikely to be foreground contaminants too.
Finally -using the CC-D for IRS 1-3 region -251 Class III-type sources including probable contaminants (14 have X-ray counterparts), 144 Class II-type sources (2 have X-ray counterparts), and 24 Class I-type YSO candidates were identified from the 'F' , 'T' , and 'P' regions, respectively. Table 4 lists these YSO candidates along with their respective NIR magnitudes, and IAU designation from Table 3 where applicable. Similar NIR CC-D could not be constructed for the IRS 9 region as we only have H and K magnitudes available for it.

Using NIR Colour-Magnitude Diagram
Many embedded and young sources can only be seen in H and K bands due to high extinction at J band wavelengths. Therefore we use the K/H-K colour-magnitude diagram (CM-D) for the identification of additional YSOs. Figure 3 shows the K/H-K CM-D for IRS 1-3, IRS 9, and the sky field regions. The almost vertical solid lines represent the zeroage-main-sequence (ZAMS) locus at a distance of 2.65 kpc reddened by AV = 0, 15, 30, 45, and 60 mag. Slanting lines indicate the reddening vectors for the marked spectral types. A low density gap in H − K colour ∼ 0.5 can be seen in all three diagrams, though it is much less pronounced for the IRS 9 region due to the lack of statistics. In the sky field CM-D (Figure 3(c)), we can see that most of the sources are confined to H − K 1. This H − K(=1) limit also corresponds to the average extinction of AV = 15 mag found towards the NGC 7538 region by Ojha et al. (2004a). Thus, in general, a background contaminant field source suffering extinction due to the cloud should also be confined to H −K 1. Also, on the same lines as Ojha et al. (2004a), if we were to assume that sources had large IR colours (say, H − K > 2) purely due to interstellar extinction, then that would imply that the reddening due to molecular cloud is AV > 30 mag. But, with such a large AV , diffuse emission will hardly be seen in NIR, which is definitely not the case here. Sources with colour over and above this limit (i.e. H − K > 1) are therefore most likely to have large H − K colour due to intrinsic infrared (IR) excess associated with YSOs. Hence, we use this colour cut-off of H − K > 1 to identify extra YSO candidates from the IRS 1-3 and IRS 9 regions. To distinguish the YSO candidates identified using the NIR CM-D from those by the CC-D (Section 3.1), we refer to these sources as 'red sources' throughout the text. 145 sources (4 X-ray counterparts) were identified in the IRS 1-3 region, and 95 (7 X-ray sources with H and K counterparts + 1 X-ray source with only K counterpart) in the IRS 9 region. The catalogs of these red sources for the IRS 1-3 and IRS 9 regions are included in Tables 4 and 5, respectively.

IRS 1-3 region
The K band image of the IRS 1-3 region is shown in Figure 4 with the overlaid YSOs and H 13 CO + (J=1-0) contours. Green squares mark the Class I sources, blue plus symbols the Class II-type or CTTS, and red circles the sources with H − K > 1. IRS 1, 2, and 3 are saturated and have been marked. Dense cloud is seen around these marked IRS sources. There appears to be relatively higher nebulosity towards the northern as opposed to the southern region. Almost all stellar sources are present in the northern region. While the Class II-type sources (blue plus symbols) seem distributed throughout the northern part, sources with H − K > 1 (red circles) and Class I sources (green squares) are mostly concentrated around the dense cloud surrounding the luminous IRS sources. Some of them could be PMS members of the embedded stellar clusters associated with the molecular clumps around IRS 1-3 and to the south of IRS 1-3. A few red sources are also seen around the ionization front (the boundary between darker nebulous region and the lighter non-nebulous region in the western part; see Section 4.3 too) at the interface between the H ii region and the molecular cloud, which might have formed due to triggered star formation. The large diffuse emission extending to the north-west of IRS sources is probably due to a combination of free-free and bound-free emissions, corresponding to what is seen optically, and coincides well with the radio brightness from the GMRT observations, while the bright and compact infrared nebula embedded within IRS 1, 2, and 3 is coincident with the peak of radio continuum (see Section 4.3). The H 13 CO + (J=1-0) contours show a peak to the south-east of the IRS 1-3 nebula, where YSO density is much lower, and most of the sources around and near this peak are the red H − K > 1 sources. Since this molecular emission traces the dense molecular cloud, it seems reasonable indeed that relatively very few sources are detected in the southern portion where the peak and the higher intensity contours lie, and the few sources detected are very reddened ones. The northern portion shows a hump in H 13 CO + (J=1-0) contours, with contours closely spaced as one moves from the northern hump towards the southern peak. Preliminary calculations were carried out to get an idea of the column densities, local thermodynamic equilibrium (LTE) mass, and virial mass of the clump associated with the peak. Assuming LTE, we tried to estimate the molecular column density using the formula from Troitsky et al. (2005), with an excitation temperature of 10 K, and a dipole moment of 3.9 debye (Botschwina et al. 1993). Fitting a 2D Gaussian to the clump results in a peak integrated line intensity of 2.0 K km s −1 and a source size at half intensity level of 1 pc. This gives us the LTE column density n(H 13 CO + ) as ∼ 3.6 × 10 12 cm −2 after applying the main beam efficiency correction. Since the abundance -X(H 13 CO + ) -estimates are variable from region to region and we do not have a measure of it here, we use the range seen for other massive star-forming regions -0.5-3.0×10 −10 from Zinchenko et al. (2009). Using these values, n(H2) range is obtained to be ∼1.2-7.2×10 22 cm −2 , and the LTE mass range of the clump to be ∼250-1500 M . The spectrum for the IRS 1-3 region clump is shown in Figure  5(upper). Applying a 1D fit to the spectrum of the clump, we obtain a linewidth of ∼ 3.4 km s −1 . The source size is corrected for the beam width to get the 'deconvolved source size' (Zinchenko 1995). Using these values along with the expression from Zinchenko et al. (1994), the virial mass is approximated to be of the order of 1000 M . Alternatively, under the reasonable assumption of gravitational equilibrium, the LTE mass and virial mass can be equated, which will imply an abundance X(H 13 CO + ) of ∼ 0.74×10 −10which is within the range derived by Zinchenko et al. (2009).
The gas-to-dust ratio for a cluster can be estimated by the relation obtained between the X-ray column density and the visual extinction for each source (Kuhn et al. 2010). Using the value of this column density (logNH1) from Table  3 and AV from IR analysis (discussed in Section 5.1.2), we find log(NH /AV ) ∼ 21.75 ± 0.20. The relation is shown in Figure 6. It is slightly higher than that for the interstellar medium, ∼ 21.34 from Ryter (1996), and other young starforming regions such as W40 (Kuhn et al. 2010). This higher value of gas-to-dust ratio is probably suggestive of dense gas in the region, as is also suggested by radio analysis (see Section 4.3).

Stellar Cluster Analysis
A stellar surface density analysis of this region was carried out using the nearest-neighbour (NN) method (Casertano & Hut 1985;Schmeja et al. 2008). We use the catalog of candidate YSOs identified in Section 3, and a procedure similar to Schmeja et al. (2008) with 20 NN. The resultant surface density map is shown in Figure 7(left), also overplotted with contours. Figure 7(right) shows the histogram of NN distances, with the peak NN distance in 0.12-0.14 range. In Figure 7(left) three major clusterings can be made out. One is towards the south, associated with the IRS 1 and IRS 3 sources, while the other two are to the north and north-west of the IRS 3 source. Each of these clusterings display multiple peaks of their own. The maximum surface density values range from 860 to 1050 pc −2 , with the maximum exhibited by the clustering to the north of the IRS 3. These values are higher than the peak calculated for the overall NGC 7538 region by Chavarría et al. (2014), which could be a result of deeper data used here. These high peak values, though, are similar to that for Serpens Core(A) (1045 pc −2 ; Schmeja et al. 2008). Figure 8 shows the K band image of the IRS 9 region with the overlaid YSOs and H 13 CO + (J=1-0) contours. Red circles denote YSO candidate sources. The main IRS sources -9, 9N1, 9N2, 9N3, and 9N4 -from Ojha et al. (2004a) have been marked on the image. As can be seen in Figure 8, sources are distributed throughout the image in no particular orientation, but with significantly more sources in the northern portion as opposed to the nebular southern portion -possibly an indication of the high extinction in this nebular region. This suggests that the region is extremely young, also affirmed by the fact that there is no free-free radio emission seen in the nebulous part (see Section 4.3). There appears to be clustered star formation going on around the IRS 9 region. The IRS 9N4 source from Ojha et al. (2004a) is resolved into two distinct sources in this image. The H 13 CO + (J=1 -0) molecular line emission peak lies to the eastern portion of the IRS 9 sources, a region deficient in YSOs. This indicates the presence of dense gas in this region. A 450 µm clump from Reid & Wilson (2005) is close to this peak. Using the same procedure as for the IRS 1-3 region (see Section 4.1), a peak integrated line intensity of 1.5 K km s −1 , and a source size at half intensity level of ∼ 0.77 pc, LTE column density n(H 13 CO + ) is calculated to be ∼ 2.7 × 10 12 cm −2 . This leads to n(H2) range of ∼ 0.9-5.4×10 22 cm −2 , and an LTE clump mass range of ∼ 100-660 M . The spectrum for the clump is shown in Figure 5(lower). The virial mass -following same steps as for IRS 1-3 region -is approximately of the order of 500 M . Alternatively, as Section 4.1, assuming gravitational equilibrium gives X(H 13 CO + ) ∼ 0.65×10 −10 .

Radio Continuum Emission
The maximum resolution radio continuum images of the entire NGC 7538 region which could be constructed are shown for 325 MHz (resolution ∼ 12.3 × 8.7 ), 610 MHz (resolution ∼ 9.0 × 4.6 ), and 1280 MHz (resolution ∼ 3 × 2 ) (Figures 9 and 10). The positions of the IRS sources in the region have been indicated by numbers in Figure 9. The radio contours in Figure 9 show a definite champagne flow morphology (Tenorio-Tagle 1979;Whitworth 1979). While the contours are closely packed in the south-west corner, they become spread as one moves from the south-west to the north-east direction. The south-west corner is density bounded, while the north-east side is ionization bounded. In the relatively lower resolution 325 and 610 MHz maps (Figure 9), there are definite peaks associated with IRS 1-3, as well as with the other IRS sources 4 and 5 (note that these IRS 4 and 5 sources are outside our NIR FoV). The peak around the IRS 1-3 region is resolved into three separate peaks -to the North, East, and West -in the higher resolution 1280 MHz image (Figure 10(left)). Out of these, the IRS 1-3-East peak coincides with the IRS 2 source. As can be seen in this high-resolution image, the corresponding radio source for IRS 2 has a cometary appearence, as has also been observed at higher radio frequencies (Bloomer et al. 1998;Campbell & Persson 1988). The IRS 1-3-North and West peaks do not have any NIR sources associated with them, though the IRS 1-3-West peak (keeping in mind that the beam size is ∼ 3 × 2 here) is very close to the IRS 2C peak from the 6 cm radio map of Campbell & Persson (1988). The IRS 1-3-East and West cores are approximately of the size of the synthesized beam (∼3 ×2 ), with a flux of ∼ 50 mJy and 30 mJy, respectively. The IRS 1-3-North core is ∼3.8 ×2.3 in size, with a total flux of ∼43.4 mJy. Based on the size of the cores (Kurtz & Franco 2002), the whole IRS 1-3 core can be classified as a compact H ii region (size from 0.1 to 0.5 pc), while the North, East, and West cores can be classified as ultracompact H ii regions (size 0.1 pc). The high resolution 1280 MHz image near the north-west part of the NGC 7538 region (near the sources marked 4 and 5 in Figure 9; outside our NIR FoV) also shows multiple cores (Figure 10(right)).
The nature of a free-free emission region can be studied by calculating its spectral index α, given by Sν ∝ ν α , where ν is the frequency and Sν is the integrated flux density of the region at ν. We can calculate the spectral index using d log Sν /d log ν. As free-free emission changes from optically thick to optically thin, α varies from 2 to -0.1 (Panagia & Felli 1975;Olnon 1975). First of all, we obtained all the GMRT images at the same uniform resolution of about ∼12 ×9 (as this is the maximum resolution map we could make for 325 MHz). Next, we calculated the integrated flux density within the IRS 1-3 compact H ii region for all three frequencies. In addition to the GMRT frequencies, VLA archival image at the frequency of 4860 MHz (at a very similar resolution of ∼14.9 ×11.5 ) was also obtained and integrated flux density calculated for the compact H ii region. The values obtained are given in Table 2. Figure  11 shows the radio SED fit using these data points. The solid grey line shows the SED fit using all four data points (α1 = 0.87 ± 0.11), while the dashed grey line shows the fitting (α2 = 1.16 ± 0.02) using only the GMRT points. These values of α suggest that the region is optically thick at these low frequencies, in consistency with earlier radio studies (Akabane et al. 1992).
If we assume that this compact H ii region is homogeneous and spherically symmetric, then we can calculate the Lyman continuum photon luminosity (photon s −1 ) using the following formula from Moran (1983,  (1) where Sν is the integrated flux density in mJy from the contour map, D is the distance in kiloparsec, Te is the electron temperature, and ν is the frequency in GHz for which the luminosity is to be calculated. Now choosing the VLA data point of 4.860 GHz (as it will lie in the most optically thin regime of the free-free emission SED), and Te a typical value of 10000 K (Panagia & Felli 1975;Olnon 1975), S * was calculated to be ∼ 1.01×10 48 (i.e. log S * ∼ 48.00). According to the tabulated values of Panagia (1973), this much flux would correspond to a ZAMS O9 star (log S * = 48.08). Using various techniques, the spectral types of IRS 1, 2, and 3 in the literature have been found to be consistent with O6 ZAMS (Willner 1976;Pestalozzi et al. 2004) (though Franco-Hernández & Rodríguez (2004 give O8.5), O9.5 ZAMS (Campbell & Thompson 1984;Bloomer et al. 1998), and a Btype star (Puga et al. 2010), respectively. However, various analyses have found extremely high extinction towards the IRS 1 source. Willner (1976) has concluded that the optical depth associated with IRS 1 is high enough to absorb 99% of UV photons. Similarly large extinction properties have been found by Campbell & Persson (1988) and Beuther et al. (2012) in their study of this source. Hence, if we exclude IRS 1, then our Lyman continuum luminosity is consistent with that for a O9.5 ZAMS source (IRS 3 contribution, being a B-type source, will be much lesser comparatively). This high extinction could be the reason why IRS 1 source is not seen in our high resolution 1280 MHz map in Figure 10(left). Additionally, this caveat of dust absorption of radiation sug-gests that O9 ZAMS spectral type should be the lower limit and the actual spectral type could be earlier than this.

LUMINOSITY AND MASS FUNCTIONS
We use the YSO catalogs from Section 3 to generate the K-band luminosity function (KLF) of the IRS 1-3 and IRS 9 stellar clusters. K band, as opposed to H or J band, is used as : the effects of extinction are minimum, it probes sources upto much fainter luminosities, and the results can be compared to the existing literature.

K−band luminosity function
KLFs for the IRS 1-3 region were generated using 2 sets of YSOs taken from Table 4 : the first set containing a combination of 'F+T+P+red-sources' (basically the entire Table 4), and a second set containing only 'T+P+Any source with X-ray detection' (see Section 3 and Figures 2 and 3). This was done as the first set may contain YSO candidates with field star contamination and will need to be corrected for it, while the second set is most likely YSOs with much lower field contamination. To derive the KLF of a region, one needs to apply corrections for the following : star count incompleteness (which is a function of magnitude), and the field star contamination towards the cluster.
Star count incompleteness was corrected for by using the completeness calculations from Section 2.1. The completeness fraction had been obtained for each 0.5 magnitude bin. The counts in each magnitude bin were scaled up by dividing the counts by the completeness fraction in the respective bins. This gives us the completeness-corrected KLF.
The field star contamination was assessed using the sky field source catalog (from Section 2.1) in conjunction with the Galactic model of Robin et al. (2003), similar to Ojha et al. (2004a,b). Two sets of catalogs were generated using the Besançon model of stellar population synthesis (Robin et al. 2003) in the direction of the sky field region. The first set contains the sources generated by setting AV = 4.4 mag, which was chosen as it is the average AV of sources in the sky field CM-D (between H − K = 0 and the low density gap at H − K = 0.5), and thus most likely to be the foreground extinction of sources towards this region. From this first catalog set, we get the total number of sources as well as the number of foreground sources in the direction of the sky field region. The foreground sources, in this catalog, are simply the ones with distances < 2.65 kpc. In addition to the total and foreground sources towards this region, we need to assess the background field star contamination too. However, the background field star contaminants will suffer an extra extinction due to the intervening molecular cloud, which needs to be taken into account. Now, the average extinction towards the NGC 7538 region has been found to be AV = 15 mag (Ojha et al. 2004a). If we assume a spherical geometry of the molecular cloud medium, then it follows that the sources behind the molecular cloud should suffer an extinction of AV = 15 + 15 = 30 mag. Therefore, we generated a second model catalog set by setting AV = 30 mag. In this second set, all the sources with distances > 2.65 kpc will give us the background contaminants. After having obtained the total, foreground, and the background sources, the K band histogram (with binwidth=0.5) was plotted and the contamination fraction (foreground+background/total) was obtained for each bin. Finally, to obtain the absolute number of contaminating sources in each magnitude bin, we use the sky field K band histogram. Each magnitude bin of the sky field K magnitude histogram is scaled by the ratio (∼ 0.75) of the areas of sky field to IRS 1-3 region, as well as by its contamination fraction calculated above. Fig. 12 shows the contamination fraction, along with the number of field contaminants in the IRS 1-3 region, as a function of magnitude. The obtained contaminant number in each bin was subtracted from the completeness-corrected KLF (binwise) to obtain the field-and completeness-corrected KLF.
The KLFs for the first YSO set ('F+T+P+redsources' ), and second YSO set ('T+P+Any source with Xray detection' ) are shown in Figure 13. A binsize of 0.5 has been chosen as it is much larger than the errors in source magnitudes. To compare with the literature, we calculate the (d log N/dmK ) slope for the rising part of the KLFs to be 0.40±0.03 (in the magnitude range 12.5-16.5) for the first YSO set, and 0.24±0.03 (in the magnitude range 12.5-16.5) for the second YSO set. While the KLF slope for the first YSO set is higher than that for the completeness-and fieldcorrected KLF for the whole NGC 7538 region (0.30 ± 0.03) from Ojha et al. (2004a), it is similar to that of the NGC 1893 (∼ 3.25 kpc; 0.34 ± 0.07) from Sharma et al. (2007), Tr 14 (∼ 2.5 kpc; 0.37 ± 0.01) from Sanchawala et al. (2007), and IRAS 06055+2039 (∼ 2.6 kpc; 0.43 ± 0.09) (Tej et al. 2006). The slope for the second set is consistent with that calculated earlier in Ojha et al. (2004a) for the completenesscorrected whole NGC 7538 region (0.28 ± 0.02) as well as that for the younger regions in NGC 7538 (0.27 ± 0.03). It is also very close to the value for the whole W3 region (∼ 1.83 kpc; 0.26 ± 0.01) (Ojha et al. 2004b) and Sh2-255 IR region (∼ 2.6 kpc; 0.17 ± 0.03) from Ojha et al. (2011). A turnoff is seen in both the KLFs after 16-16.5 magnitude bin, similar to the Tr 14 region from Sanchawala et al. (2007).

Mass Function
The mass function (MF) is usually described by the following differential form : where N is the number of stars and M * is mass of stars (Bastian et al. 2010;D'Antona 1998;Chabrier 2003;Kroupa et al. 2013). Since magnitude (mK ), and not mass (M * ), is an observable quantity, ξ has to be re-written as : The absolute magnitudes were converted to apparent magnitudes using the distance modulus (at a distance of 2.65 kpc). The complete MLR is shown in Figure 14, along with the curve fit to it. To calculate ξ, we need the derivative of this MLR. We limit our analysis to sources with mK 12 mag, as it the lower limit of the MLR.
In addition to the MLR, a form of the intrinsic KLF is also needed. To derive the intrinsic KLF, the sources were corrected for extinction by dereddening them along the reddening vectors to the respective loci as per the following order of steps. First, the sources which were present in the 'F' region of the NIR CC-D (Figure 2(a)) were dereddened to the dwarf locus (whose low-mass regime -from the turnover onwards -was approximated by a straight line; similar to Tej et al. 2006;Samal et al. 2007), while those in the 'T' and 'P' regions were dereddened to the CTTS locus. The rest of the sources (those not dereddened using the NIR CC-D) were dereddened using the H and K magnitudes to the ZAMS locus (approximated by a straight line) in the NIR CM-D (Figure 3(a)) to estimate their extinctions, though this is only an approximation as the YSOs are most-likely not on ZAMS yet. The major source of uncertainty in using the CM-D is the correction of the locus for distance, which need not be well-determined. The individual visual extinction of the sources was found to range upto 60 mag. This is much deeper, as expected, than upto 40 mag which was estimated by Ojha et al. (2004a) for their detected YSOs. The histogram (not shown here) of these visual extinctions was found to peak in the range 7.5-10 mag. The mean extinction value was ∼ 15 mag, consistent with that used for model simulations in Section 5.1.1 and Ojha et al. (2004a). The sources with AV > 40 mag were mostly found to be distributed along the southern boundary between the nebulous and non-nebulous region (also see Figs. 1(left) and 4). A few such sources were also associated with the small nebular patch at the south-east corner of the IRS 1-3 region.
Using this catalog of reddening-corrected sources, field-, completeness-, and reddening-corrected KLFs were obtained for the two YSO sets ('F+T+P+red-sources' being the first set and 'T+P+Any source with X-ray detection' being the second set; see Section 5.1.1). To preserve the information about the form of the KLF, different magnitude intervals (each interval containing three or more bins) were fit with equations of straight line using simple linear regression. This gives us the KLF slope for each magnitude bin.
After we obtain the KLF slope (dN/dmK ) in each magnitude bin and the form of the MLR derivative (d log M * /dmK at discrete points), value of ξ is calculated using Equation 3 for each magntiude bin -which is further mapped onto the log M * space using the MLR. The resulting form of ξ is shown in Figures 15 and 16 along with the KLFs. Poissonian error of ± √ N is marked for each bin. It should be kept in mind that equal magnitude intervals do not map onto equal log M * intervals. The shape of the ξ is closest to that derived by Scalo (1986) for field stars, with a peak at the low mass end, and another peak at the intermediate mass (Meyer et al. 2000). The low mass peak is at log M * bin of −1.09 to −1.00, i.e 0.08 − 0.1 M (for Figures  15 and 16). This rise in the MF till the BD limit has been seen for other regions like W3 Main (Ojha et al. 2009) and S106 (Oasa et al. 2006) too, as well as mentioned likely by Kroupa (2007). The curve does not extend enough in the intermediate mass range to clearly discern the peak.
Another way to test the MF form is to use the MLR to assign a mass to each star in the catalog, and then bin those masses to obtain the N (log M * ), rather than the differential form ξ(log M * ) (Chabrier 2003, see its Section 1.3). Using this, we obtain the MFs shown in Figure 17. The field star subtraction for the first catalog here was carried out statistically using the NIR CM-D as follows. Similar to Sharma et al. (2007), we divided the CM-Ds of the IRS 1-3 region and the sky field region into grids with ∆K=0.5 mag and ∆(H − K)=0.1 mag. The number of stars was compared on a grid-by-grid basis, and from the IRS 1-3 CM-D, a fixed number of stars (equal to the number in the corresponding sky field grid) was removed based on their distances to the sky field stars in the colour-magnitude space. Each grid was corrected for incompleteness.
As can be seen, the MFs are similar in shape to those in Figures 15 and 16  . A secondary peak seen here is also observed in the IMF of the ρ Ophiuchi cluster from Alves de Oliveira et al. (2012) and Orion Nebula Cluster from . If we assume a power law form of the MF (similar to Salpeter 1955), then : Now, if we additionally assume that the star formation is strictly coeval, then it can be mathematically shown that and that the present day mass function will have the same slope as the IMF (Massey 1998, see its Section 2.1). Also, since the age of this region is ∼ 1 Myr, which is lower than the MS life of even the most massive stars, none of the stars would have disappeared from the field. In Figure 17, fitting the sub-solar low mass range (between the first and the second peaks, i.e. 0.1-1 M ) using Equation 6, we get the value of Γ for the first and second YSO sets as 0.54 ± 0.05 (say, Γ1) and 0.75 ± 0.18 (say, Γ2), respectively. Both the slopes are lower than the Salpeter slope of 1.35. While Γ1 seems consistent with that from Kroupa (2002) (also see Figure 2 of review by Bastian et al. 2010), Γ2 is slightly steeper. This steepness, though, might be explained by the fact that Figure 17(right) includes mostly the youngest sources with fainter magnitudes -thus making this MF 'bottom-heavy' , i.e. more sources at lower mass ranges. On the other hand, Figure 17(left) includes Class III sources and thus more sources in relatively higher mass ranges. In general, more complicated/realistic forms of the IMF can emerge due to accretion processes, leading to a tail towards the high mass end (Dib et al. 2010), which will become noticeable only in high-mass star-forming regions where the intermediate to high-mass bins are well populated. Finally, we should keep in mind the caveat that any derivation of an IMF suffers from multiple unavoidable and systematic biases, e.g. those arising due to -among others -choice of the stellar MLR, PMS evolution, effect of unresolved sources, binning, etc. Kroupa et al. (2013, see their Section 2.1) have dealt with these biases in a succinct manner. Different treatment of these biases might accordingly alter a derived IMF. The ratio of stars to BDs in a region is often used as another quantitative indicator for the mass function. We consider the following two definitions of this ratio :

IRS 9 region
We obtained the KLF for the IRS 9 region using the catalog of sources identified in at least K band, as well as the catalog of candidate YSOs (red sources with H − K > 1 identified in Section 3.2 and those with X-ray emission). The KLFs are shown in Figure 18. Figure 18(left) shows the raw KLF (grey line) along with the field-and completeness-corrected KLF (black line). Since the statistics are low for this region, we only make a few qualitative comparisons here. The (d log N/dmK ) slope was calculated to be 0.20 ± 0.04 (in 13-19 K mag range), and is lower than that calculated for younger regions in Ojha et al. (2004a) and some other starforming regions (see Section 5.1.1), though similar to Sh2-255 IR region (Ojha et al. 2011). The KLF of red sources (H − K > 1) is shown in Figure 18(right). There appear to be two peaks in Figure 18(right) -one near 15-15.5 mag and another at about the completeness limit (marked by a dotted vertical line). The KLF seems to be plateauing near the completeness limit. The brighter peak is most likely to be due to field contamination, which has not been corrected for due to a lack of statistics here. If we consider the 16-17.7 mag interval of this YSOs' KLF, we see that it first has a steeper slope in the 16-16.5 mag interval, and then plateaus off in the 16.5-17.7 mag interval. This steep slope followed by plateauing is seen in other regions like NGC 1624 (Jose et al. 2011) and Sh2-255 IR region (Ojha et al. 2011) too, albeit for different magnitude limits. In this range (16-17.5), the slope comes out to be 0.24 ± 0.11, which is very similar to that from Ojha et al. (2004a) for younger regions of NGC 7538, and is consistent within errors with that obtained for the Sh2-255 IR region (Ojha et al. 2011).

DISCUSSION
Analysis of star-forming regions in different stages is essential in understanding how star formation proceeds and the effect of manifold physical processes on various diagnostic tools (like LF, MF, etc). In general, however, physical conditions (like density, temperature, chemical composition) could differ from one region to another. The study of NGC 7538 cluster regions (IRS 4-6, IRS 1-3, and IRS 9, which follow an age sequence in descending order; Ojha et al. 2004a), being part of the same complex, should mitigate this problem. Here, using the deepest NIR data, we have obtained the KLF and the MF of the IRS 1-3 region, while only the KLF of the IRS 9 region has been discussed. The MF for the IRS 1-3 region shows that it rises till the BD limit before turning over, which indicates lower temperature or denser gas distribution than for Orion nebula cluster Muench et al. 2002), though comparison of LF and MF of different regions from literature always follow the caveat that the method used for their construction might be slightly different. Presence of dense molecular material has been attested by surveys in submillimetre ranges too (Sandell & Sievers 2004;Reid & Wilson 2005). Our molecular line analysis also reveals massive clumps which could fragment and lead to future stellar cluster formation. The 850 µm emission (Chavarría et al. 2014) shows the IRS 1-3 region to be located on a junction of filaments, which could partially explain the active star formation going on. Confinement of most IR clusters at filament junctions has also been observed by Schneider et al. (2012) and obtained in the simulations of Dale, Ercolano, & Bonnell (2012). The entire morphology of NGC 7538 seems similar to the hub-filament structure proposed by Myers (2009), where the IRS 1-3 cluster region forms the 'hub' to which various filaments merge. The molecular hydrogen column density (which is > 10 23 cm −3 ; Fallscheer et al. 2013) satisfies the condition which Myers (2009) proposes for a 'hub' . The IRS 9 region, however, appears to be distinct and at one end of the filament which connects it with IRS 1-3 (this filament goes on to connect IRS 4 too; see Sandell & Sievers 2004). IRS 9 region has been found to be much younger than IRS 1-3, owing to a lack of free-free radio emission as well as a low number of YSOs associated with the nebula. The Herschel temperature maps (Fallscheer et al. 2013) also show that IRS 9 is much colder than IRS 1-3 region. The decreasing age-sequence along the filament connecting the three cluster regions of IRS 4-6, IRS 1-3, and IRS 9 alludes to the possibility that star formation could have been triggered along this filamentary structure.
To gain a firm understanding of the NGC 7538 region stellar population, as well as the evolution of the KLF and MF as one moves from relatively older to relatively younger regions in a star-forming region, we need to obtain the detailed LFs and MFs for the IRS 4-6 region as well as IRS 9 region. The advantage will be that all the sub-regions belong to the same region NGC 7538, and should therefore differences arising due to the differences in physical conditions (leading to different temperature and density distributions, which leads to different Jeans masses and thus different LFs/MFs) need not be a concern. Future analysis of the BDs detected in the region, construction of MF from a spectroscopic sample, an analysis of core mass function of this region, and examination of spectral energy distributions of individual sources are also required for a better idea of the star formation going on.

CONCLUSIONS
We have carried out the deep NIR imaging survey of IRS 1-3 (J, H, and K) and IRS 9 (H and K) sub-regions in the NGC 7538 star-forming region with the highest spatial resolution so far. In addition, GMRT observations at 325, 610, and 1280 MHz are used to examine the radio emission and physical characteristics. H 13 CO + (J=1-0) molecular line emission from Nobeyama radio telescope is used to understand the morphology. Our main results are summarized as follows.
(i) Based on the NIR CC-D (J − H/H − K), 144 Class IItype, and 24 Class I-type YSOs were identified in the IRS 1-3 region. Using the NIR CM-D (K/H − K), 145 sources were identified in the IRS 1-3 region, and 96 in IRS 9 region. 27 sources were found to have X-ray counterparts.
(ii) In the IRS 1-3 region, the red sources and the Class I sources are concentrated around the compact nebula associated with luminous IRS sources, and the H 13 CO + (J=1-0) contours show a peak (column density n(H 13 CO + ) ∼ 3.6 × 10 12 cm −2 ) to the south of these sources. Stellar surface density analysis reveals three clusterings in this region. The IRS 9 region does not have any particular distribution of the YSOs, with the nebula hardly containing any sources. An H 13 CO + (J=1-0) peak (n(H 13 CO + ) ∼ 2.7 × 10 12 cm −2 ) lies to the east of the cluster around IRS 9 source. The virial masses are approximately of the order of 1000 M and 500 M for the clumps in IRS 1-3 and IRS 9 regions, respectively.
(iii) Radio emission shows a champagne flow morphology in the NGC 7538 region. In low-resolution 325 and 610 MHz maps, a compact H ii region is seen associated with IRS 1-3 sources, with its spectral index calculated to be 0.87 ± 0.11, suggesting optically thickness. In very high resolution 1280 MHz maps, the IRS 1-3 compact H ii region is resolved into three separate peaks, one of which coincides with the known IRS 2 source.
(iv) KLFs were constructed for IRS 1-3 region for two sets of sources : 'F+T+P+red-sources' and 'T+P+Any source with X-ray detection' . The rising part of the KLF has a (d log N/dmK ) slope of 0.40±0.03 for the first set, and 0.24± 0.03 for the second set (K magnitude range 12.5-16.5). The KLF slopes for the IRS 9 region using sources with K band only detections was calculated to be 0.23±0.04 (13-19 mag).
(v) Theoretical mass-luminosity relation is used to obtain ξ(log M * ) (differential form of MF) and N (log M * ) for the IRS 1-3 cluster region. Both ξ(log M * ) and N (log M * ) show a peak in the low mass regime as well as a peak in the intermediate mass regime. In low-mass regime, ξ(log M * ) extends upto BD limit before turn-off (0.08-0.1 M ), while N (log M * ) peak is in 0.1-0.18 M range. The slope d log N/d log M * (from N (log M * ) distribution) for the first and second YSO sets in the range 0.1-1 M are 0.54 ± 0.05 and 0.75 ± 0.18, respectively -much lower than the Salpeter value of 1.35. The MFs most closely resemble that of Scalo (1986). The star to BD ratio upper limit was calculated to be 10.2.
(3)  The IAU designation is from Table 3. Table 4 is available in its entirety in a machine-readable form in the online journal. A portion is shown here for guidance regarding its form and content. The IAU designation is from Table 3. Table 5 is available in its entirety in a machine-readable form in the online journal. A portion is shown here for guidance regarding its form and content.   Bessell & Brett (1988). The CTTS locus (Meyer et al. 1997) is shown by the blue solid line, and its extension into the 'P' region by grey dot-dashed line. The three parallel and slanted dashed lines are the reddening vectors using the reddening laws of Cohen et al. (1981). Blue crosses mark the sources with X-ray counterparts. All points are in CIT system. (b) The CC-D for the sky field. Dec (2000) 3 2 1 Figure 4. K band image of the IRS 1-3 region overlaid with YSOs and H 13 CO + (J=1-0) contours (black lines). The contours are drawn in the range 0.5-2.2 K km s −1 , with step size of 0.1 K km s −1 . Class I sources are shown with green square symbols, Class II with blue plus symbols, and sources with H − K > 1 with red circles. IRS 1, 2, and 3 are saturated in the image and have been marked. Figure 5. The H 13 CO + (J=1-0) molecular line spectra towards (upper) the peak in IRS 1-3 region, and (lower) the peak in IRS 9 region.  Figure 9. The maximum resolution radio continuum images for 325 MHz (cyan contours; 12.3 × 8.7 ) and 610 MHz (green contours; 9.0 × 4.6 ) obtained using GMRT for the entire NGC 7538 region. The contours for 325 MHz image are drawn at 12,17,25,30,35,37,40, and 42 σ (σ ∼ 1.78 mJy). 610 MHz contours are drawn at 2.5, 5, 7, 10, 15, and 20 σ (σ ∼ 4 mJy). The IRS sources in this FoV are marked by their respective numbers in magenta. The background image is 325 MHz radio image.    Robin et al. (2003). (lower) The number of contaminants (per magnitude bin) for the IRS 1-3 cluster region. The error bars show the error due to counting statistics. A similar exercise was carried out for the IRS 9 region.